Skip to main content

Fused ultrasound and electromyography-driven neuromuscular model to improve plantarflexion moment prediction across walking speeds

Abstract

Background

Improving the prediction ability of a human-machine interface (HMI) is critical to accomplish a bio-inspired or model-based control strategy for rehabilitation interventions, which are of increased interest to assist limb function post neurological injuries. A fundamental role of the HMI is to accurately predict human intent by mapping signals from a mechanical sensor or surface electromyography (sEMG) sensor. These sensors are limited to measuring the resulting limb force or movement or the neural signal evoking the force. As the intermediate mapping in the HMI also depends on muscle contractility, a motivation exists to include architectural features of the muscle as surrogates of dynamic muscle movement, thus further improving the HMI’s prediction accuracy.

Objective

The purpose of this study is to investigate a non-invasive sEMG and ultrasound (US) imaging-driven Hill-type neuromuscular model (HNM) for net ankle joint plantarflexion moment prediction. We hypothesize that the fusion of signals from sEMG and US imaging results in a more accurate net plantarflexion moment prediction than sole sEMG or US imaging.

Methods

Ten young non-disabled participants walked on a treadmill at speeds of 0.50, 0.75, 1.00, 1.25, and 1.50 m/s. The proposed HNM consists of two muscle-tendon units. The muscle activation for each unit was calculated as a weighted summation of the normalized sEMG signal and normalized muscle thickness signal from US imaging. The HNM calibration was performed under both single-speed mode and inter-speed mode, and then the calibrated HNM was validated across all walking speeds.

Results

On average, the normalized moment prediction root mean square error was reduced by 14.58 % (\(p=0.012\)) and 36.79 % (\(p<0.001\)) with the proposed HNM when compared to sEMG-driven and US imaging-driven HNMs, respectively. Also, the calibrated models with data from the inter-speed mode were more robust than those from single-speed modes for the moment prediction.

Conclusions

The proposed sEMG-US imaging-driven HNM can significantly improve the net plantarflexion moment prediction accuracy across multiple walking speeds. The findings imply that the proposed HNM can be potentially used in bio-inspired control strategies for rehabilitative devices due to its superior prediction.

Background

Locomotion mobility accounts for a dominant part of human activities of daily living, like moving around the home and community, going to work or school, doing errands, visiting friends, etc. The human lower extremity plays an essential role in achieving locomotion mobility. The human ankle plantarflexors generate a large burst of “push-off” mechanical power during the late stance phase of walking, enabling forward and upward acceleration of the body’s center of mass. Due to neurological disorders or injuries like spinal cord injury, stroke, and multiple sclerosis, the weakened function or dysfunction of plantarflexors is likely to cause a dramatic decrease in the “push-off” power. Consequently, these mobility disorders impair walking function and cause poor energy economy [1], as well as disrupt both physical and emotional well-being [2].

Modern neurorehabilitation devices, such as powered ankle exoskeletons [3,4,5], soft exosuits [6, 7], and functional electrical stimulation [8,9,10,11,12], may use assist-as-needed control to actively engage and maximize recovery of users with mobility impairments [13, 14]. In turn, the efficacy of the control strategy depends on the accurate determination of continuous human volitional movement intent (net joint moment). Mechanical sensors, like force or torque sensors, installed on a rigid and bulky frame, have been often used to measure the human intent for joints without direct interaction to the ground, but limit the system’s wearability. Also, inaccuracies may creep in easily due to the inevitable misalignment between the bionic joint center and human joint center, which may introduce undesired interaction force [15, 16]. Usually, it is challenging to directly measure the net ankle joint moment during walking overground with conventional force or torque sensors setup. The standard way to measure the net ankle joint moment uses a motion capture system, ground reaction force (GRF), and inverse dynamics (ID) calculations. However, there are two shortcomings of the standard approach. First, the setup is constrained to a lab environment, and not applicable for field testing. Second, the results from ID do not reveal how skeletal muscles perform during walking from a neuromuscular perspective. Therefore, a forward dynamics approach based on a neuromuscular model would be very instrumental when a motion capture system and GRF data are unavailable.

Surface electromyography (sEMG) measures electrical potentials during asynchronous muscle neurons firings, and its amplitude and frequency positively relate to muscle activation levels. Therefore, sEMG-derived signals can be used in a Hill-type neuromuscular model (HNM) or to train a machine learning approach (model-free) to predict volitional joint moment [17, 18] and angular position [19,20,21]. However, sEMG signals suffer from interference or cross-talking from the adjacent muscles, and an inability of measuring activations of deep-layer muscles [22, 23]. Alternatively, two-dimensional brightness mode (B-mode) ultrasound (US) imaging allows one to see the musculature of the targeted muscle in vivo. Due to its ability to directly visualize superficial and deep-layer muscles, US imaging may work as an alternative methodology to predict joint motion or motion intent. Potentially, US imaging overcomes the shortcomings of the sEMG measurements. Most frequently used architectural features from US images include pennation angle (PA) [24, 25], fascicle length (FL) [26, 27], muscle thickness (MT) [28, 29], and cross-sectional area [30]. These features have been correlated with the joint kinetics and kinematics during isometric or isokinetic joint motion by using HNM-based or model-free approaches [25, 31,32,33].

Motivation also exists to use a dual-modal approach that combines the measurement from sEMG and US imaging. Potentially sEMG signals and US imaging provide complementary information, and the combination between them may (1) mitigate any cross-talking effect from neighboring sEMG signals and (2) lower US imaging-derived features’ drift due to the accumulated errors from cyclic joint movement. Our recent studies have shown the advantages of using the dual-modal approach over uni-modal bio-signals (sEMG or US imaging) for ankle joint moment/motion prediction under isometric/dynamic dorsiflexion studies [21, 25, 27] and isometric plantarflexion study [34, 35]. Similarly, Dick et al. used both sEMG and US imaging to predict plantarflexion force during dynamic cycling tasks [32]. However, the use of the dual-modal bio-signals during more complex functional tasks, like walking across different speeds remains unexplored.

In this study, we investigated a dual-modal approach that takes both processed sEMG and MT from US imaging as inputs to a modified HNM, named sEMG-US imaging-driven HNM, to predict net ankle joint plantarflexion moment during the walking stance phase across multiple speeds. We hypothesize that (1) the proposed HNM will achieve a better net plantarflexion moment prediction accuracy than sEMG-driven and US imaging-driven HNMs, (2) the net plantarflexion moment prediction performance is more robust if data collected from multiple speeds are included in HNMs’ calibration procedure.

In previous US imaging studies, visualized architectural features, such as PA and FL, require high-resolution US imaging, which can be significantly affected by a US transducer placement site on the muscle. To mitigate the requirement of high-resolution US imaging, US imaging-derived signals such as echogenicity [27, 36, 37], tissue displacement [38], and MT [29], are more preferable to correlate with muscle or joint mechanical functions. Because ankle plantarflexors: lateral and medial gastrocnemius and soleus muscles (LGS, MGS, and SOL) are not accessible in the same US image plane, this study chose to focus on LGS and SOL, which are in the same plane, and tracked their MT change during the walking experiments across multiple speeds. Therefore, the contributions of the paper are: (1) the use of US imaging-derived MT changes and sEMG-derived changes of LGS and SOL muscles as surrogate variables of muscle activations during walking, (2) development of a modified HNM that uses both sEMG and US imaging as inputs for single-speed modes and an inter-speed mode, and (3) evaluation of the sEMG-US imaging-driven HNM’s robustness across multiple walking speeds.

Method

Hill-type neuromuscular model of ankle joint

Below, we propose a modified sEMG-US imaging-driven HNM to directly build a relationship between the joint net PF moment and sEMG-US imaging-derived surrogate signals of both LGS and SOL muscles. There are three sub-models involved in the HNM: (a) sEMG-US imaging-derived weighted muscle activation model, (b) muscle-tendon unit geometry model, and (c) muscle contraction dynamic model.

Weighted muscle activation model

The neural activation at \(t_{k}\) time instant \(N_{i}(t_{k}),\,i=1,2\), \(k=1,2,3,...\), for SOL and LGS muscles, respectively, considers the electromechanical delay (EMD), \(\tau\), between the onset of an sEMG signal and a muscle contraction and utilizes a second-order recursive filter and is defined as [39]

$$\begin{aligned} N_{i}(t_{k})=\alpha _{i}u_{i}(t_{k-\tau })-\beta _{1i}N_{i}(t_{k-1})-\beta _{2i}N_{i}(t_{k-2}) \end{aligned}$$
(1)

where EMD, \(\tau\), is usually between 30 ms and 120 ms. \(u_{i}(t_{k})\) represents the sEMG’s linear envelope normalized to the peak value of the specific task (with a constant peak value across all walking speeds on each subject in this study). The sEMG’s linear envelope was derived after raw sEMG signals’ band-pass filtering, full-wave rectification, and low-pass filtering. \(\alpha _{i},\) \(\beta _{1i}\), and \(\beta _{2i}\) are coefficients that define the recursive filter’s dynamics of each muscle [39], and the following set of constrains is employed to reach a positive stable solution, i.e.

$$\begin{aligned} \beta _{1i}=\gamma _{1i}+\gamma _{2i},\,\beta _{2i}=\gamma _{1i}\cdot \gamma _{2i},\,\alpha _{i}-\beta _{1i}-\beta _{2i}=1 \end{aligned}$$
(2)

where \(\left| \gamma _{1i}\right| <1\) and \(\left| \gamma _{2i}\right| <1\).

A nonlinear relationship between neural activation \(N_{i}(t_{k})\) and corresponding muscle activation \(a_{1i}(t_{k})\) is given as [39]

$$\begin{aligned} a_{1i}=\frac{e^{A_{i}N_{i}}-1}{e^{A_{i}}-1} \end{aligned}$$
(3)

where \(A_{i}\) represents the nonlinear shape factor for each muscle, which is allowed to vary between − 3 and 0, with \(A_{i}=-3\) being a nonlinear and \(A_{i}=0\) being a linear relationship.

According to the reported results in [40,41,42,43], the MT change was found to correlate with the muscle contraction level or muscle activation through a linear function or piece-wise linear function. Therefore, in this work, the second part of the LGS or SOL muscle activation is calculated from the US imaging-derived MT change. The MT values of LGS and SOL muscles are denoted as \(MT_{i}(t_{k})\). According to the preliminary results in [44], there is a positive relationship between MT change and the targeted muscle contraction level during the walking stance phase. After taking normalization of the MT with respect to the lower bound (at rest state) and upper bound (peak value of the specific task) on each participant, the US imaging-derived muscle activation, \(a_{2i}(t_{k})\), is proposed as

$$\begin{aligned} a_{2i}=\frac{MT_{i}-MT_{i\min }}{MT_{i\max }-MT_{i\min }} \end{aligned}$$
(4)

where \(MT_{i\max }\) and \(MT_{i\min }\) represent the constant subjective MT values when the LGS and SOL muscles are at the walking task-specific maximum voluntary contraction and complete rest condition, respectively. For each individual, \(MT_{i\max }\) and \(MT_{i\min }\) are set as consistent values across different walking speeds. This normalization guarantees that the US imaging-derived muscle activations vary between 0 and 1.

By introducing an allocation gain between the sEMG- and US imaging-derived muscle activations, the synthesized/weighted muscle activation levels for LGS or SOL muscle, \(a_{i}(t_{k})\), can be represented as

$$\begin{aligned} a_{i}=\delta _{i}a_{1i}+(1-\delta _{i})a_{2i} \end{aligned}$$
(5)

where \(\delta _{i}\in [0,\,1]\) represents the muscle activation allocation gain for LGS (\(i=1\)) and SOL (\(i=2\)) muscles.

Muscle-tendon unit geometry model

Fig. 1
figure 1

Schematic diagram of the muscle-tendon-unit (MTU) geometry model and contraction dynamic model of the plantarflexor muscles during treadmill walking tasks. a A human participant walking on a treadmill. b The MTU geometry model, where the ankle joint rotation center in the sagittal plane is denoted as O, and the proximal and distal osteotendinous junction points near the knee joint and heel are denoted as A and B. SOL and LGS muscles are referred as \(i=1,2\). c The MTU contraction dynamic model, where the anatomical structure and contraction force generation of the MTU are illustrated

As presented in Fig. 1b, consider the ankle joint rotation center in the sagittal plane as O, the proximal and distal osteotendinous junction points of each each muscle-tendon unit (MTU) near the knee joint and heel as \(A_{i}\) and \(B_{i}\), and the angle between \(OA_{i}\) and \(OB_{i}\) as \(q(t_{k})\), then each MTU length, \(l_{mt_{i}}(t_{k})\), is represented as the distance between \(A_{i}\) and \(B_{i}\), and is calculated as

$$\begin{aligned} l_{mt_{i}}=\sqrt{l_{OA_{i}}^{2}+l_{OB_{i}}^{2}-2l_{OA_{i}}l_{OB_{i}}\cos (q)} \end{aligned}$$
(6)

where \(l_{OA_{i}}\) and \(l_{OB_{i}}\) represent the distance of \(OA_{i}\) and \(OB_{i}\) obtained from OpenSim (National Institutes of Health for Biomedical Computation, Stanford, USA) [45], respectively. A generic OpenSim model (gait2392) was linearly scaled to each participant in OpenSim version 4.1 [45] per [46]. According to the law of sines, the moment arm of each MTU, \(r_{mt_{i}}(t_{k})\), is calculated as

$$\begin{aligned} r_{mt_{i}}=\frac{\partial l_{mt_{i}}(q)}{\partial (q)}=\frac{2l_{OA_{i}}l_{OB_{i}}\sin (q)}{l_{mt_{i}}}. \end{aligned}$$
(7)

The MTU generates contraction force only when it is stretched, which indicates the current tendon length \(l_{t_{i}}(t_{k})\) is equal to or longer than the tendon slack length \(l_{t_{i}}^{sk}\). From the perspective of muscle contraction dynamics that is shown in Fig. 1c, the overall MTU length, \(l_{mt_{i}}(t_{k})\), could also be expressed as

$$\begin{aligned} l_{mt_{i}}=l_{t_{i}}+l_{m_{i}}\cos (\phi _{i}) \end{aligned}$$
(8)

where \(l_{t_{i}}(t_{k})\) and \(l_{m_{i}}(t_{k})\) represent the current tendon length and muscle fascicle length for LGS and SOL muscles, respectively. \(\phi _{i}(t_{k})\) represents the pennation angle that changes with instantaneous \(l_{m_{i}}(t_{k})\). By assuming the individual muscle belly has a relatively small MT change and volume change [47, 48], \(\phi _{i}(t_{k})\) can be approximately calculated as

$$\begin{aligned} \phi _{i}\approx \arcsin \left( \frac{l_{m_{i}}^{0}\sin \phi _{i}^{0}}{l_{m_{i}}}\right) \end{aligned}$$
(9)

where the term \(\phi _{i}^{0}\) denotes the constant pennation angle when the muscle is at optimal fascicle length \(l_{m_{i}}^{0}\). The work in [49, 50] has shown that the optimal fascicle length increases as muscle activation decreases. Due to this coupling, the following relationship between the muscle activation and corresponding optimal fascicle length is given as

$$\begin{aligned} l_{m_{i}}^{0}=l_{mo_{i}}^{0}(\lambda (1-a_{i})+1) \end{aligned}$$
(10)

where \(\lambda\) is the rate of change in the optimal fascicle length, and it is selected as 0.15 [39]. \(l_{mo_{i}}^{0}\) is the optimal fascicle length at the walking task-specific maximum voluntary contraction, and \(l_{m_{i}}^{0}(t_{k})\) is the optimal fascicle length at time \(t_{k}\) and muscle activation \(a_{i}(t_{k})\).

From the preliminary results of plantarflexors’ US imaging in [44], it is very challenging to capture the entire fascicle length \(l_{m_{i}}(t_{k})\) of LGS or SOL muscles by using the current US transducer with a small dimension (width of 38 mm). Therefore, the fascicle length parameters of either LGS or SOL muscles were not measured directly from the US images. Instead, an indirect calculation method was applied as detailed below. According to [51], the nonlinear nominal tendon force-tendon strain relationship is given as

$$\begin{aligned} F_{t_{i}}(\xi _{i})=\left\{ \begin{array}{ll} 0, &{} \xi \le 0\\ 1480.3F_{i}^{\max }\xi _{i}^{2}, &{} 0<\xi <0.0127\\ (37.5\xi _{i}-0.2375)F_{i}^{\max }, &{} \xi \ge 0.0127 \end{array}\right. \end{aligned}$$
(11)

where \(\xi _{i}=\frac{l_{t_{i}}(t_{k})-l_{t_{i}}^{sk}}{l_{t_{i}}^{sk}}\) represents the tendon strain and \(F_{i}^{\max }\) represents the muscle contraction force at the walking task-specific maximum voluntary contraction. By considering \(F_{t_{i}}(\xi _{i})=F_{mt_{i}}\), where \(F_{mt_{i}}\) is defined in (13), \(l_{t_{i}}(t_{k})\) can be numerically computed based on the Runga–Kutta–Fehlberg algorithm. By substituting (6) and (9) to (8), the muscle fascicle length \(l_{m_{i}}(t_{k})\) will be calculated, as well as the muscle fascicle velocity \(v_{m_{i}}(t_{k})\) by taking the time derivative of \(l_{m_{i}}(t_{k})\).

Muscle dynamic contraction model

In the sEMG-US imaging-driven HNM, the individual moment component produced by each muscle, \(M_{i}(t_{k})\), can be represented as

$$\begin{aligned} M_{i}=F_{mt_{i}}r_{mt_{i}} \end{aligned}$$
(12)

where \(r_{mt_{i}}(t_{k})\) is defined in the above geometry model and \(F_{mt_{i}}(t_{k})\) denotes the corresponding contraction force generated on each MTU and is represented as

$$\begin{aligned} F_{mt_{i}}=(F_{ce_{i}}+F_{pe_{i}})\cos (\phi _{i}) \end{aligned}$$
(13)

where \(\phi _{i}(t_{k})\) represents the pennation angle defined above. \(F_{ce_{i}}(t_{k})\) and \(F_{pe_{i}}(t_{k})\) denote the corresponded forces generated by the parallelly located contractile element and the passive element, and can be calculated as

$$\begin{aligned} \begin{array}{lll} F_{ce_{i}} &{} = &{} F_{i}^{\max }f_{l_{i}}(l_{m_{i}})f_{v_{i}}(v_{m_{i}})a_{i}\\ F_{pe_{i}} &{} = &{} F_{i}^{\max }f_{p_{i}}(l_{m_{i}}) \end{array}. \end{aligned}$$
(14)

In (14), \(F_{i}^{\max }\) of the LGS and SOL muscles will be identified based on optimization algorithm in the HNM calibration procedures. \(f_{l_{i}}(l_{m_{i}}(t_{k}))\), \(f_{v_{i}}(v_{m_{i}}(t_{k}))\), and \(f_{p_{i}}(l_{m_{i}}(t_{k}))\) denote the generic muscle contractile force-fascicle length, force-fascicle velocity, and passive elastic force-fascicle velocity curves. These curves were normalized to \(F_{i}^{\max }\), optimal fascicle length \(l_{m_{i}}^{0}\), and maximum fascicle contraction velocity \(v_{m_{i}}^{\max }\). The explicit expressions of \(f_{l_{i}}(l_{m_{i}}(t_{k}))\), \(f_{v_{i}}(v_{m_{i}}(t_{k}))\), and \(f_{p_{i}}(l_{m_{i}}(t_{k}))\) can be found in [25, 32, 52].

Experimental protocol, data collection and pre-processing

Fig. 2
figure 2

The treadmill walking experimental setup, data collection, processing, and Hill-type neuromuscular model calibration. a Illustration of walking experimental setup. (1) Instrumented treadmill with two split belts and in-ground force plates. (2) 39 retro-reflective markers on participant’s lower body for lower limb kinematics measurements. (3) Four sEMG sensors to record signals from LGS, MGS, SOL, and tibialis anterior (TA) muscles. (4) A single US transducer to image both LGS and SOL muscles with the appropriate probe placement. (5) Ultrasound imaging machine to collect radio frequency data. (6) Computer screen to show B-mode ultrasound imaging. (7) Computer screen to show live markers and segment links of the participant. (8) 12 motion capture cameras to track markers’ trajectories. b Schematic illustration of the proposed sEMG-US imaging-driven HNM calibration. The weighted muscle activation for LGS or SOL combines both sEMG signals and US imaging. Red solid, black solid, and blue dashed lines with an arrow represent the input signals to HNMs, the intermediate sub-model output signals, and the unknown parameters adjustment, respectively

The study was approved by the Institutional Review Board (IRB) at North Carolina State University (IRB number: 20602). Ten young participants (7M/3F, age: 25.4 ± 3.1 years, height: 1.77 ± 0.10 m, mass: 78.0 ± 21.1 kg) without any neuromuscular or orthopedic disorders, were recruited in this study. Every participant got familiar with the experimental procedures and signed an informed consent form before participation.

Figure 2a summarizes the experimental setup for the participants to perform static anatomical poses and dynamic gait trials (walking speeds changing from 0.50 to 1.50 m/s), and Fig. 2b presents the workflow of the data processing and sEMG-US imaging-driven HNM calibration procedures. During all walking trials, three-dimensional coordinates of 39 retro-reflective markers positioned on the participant’s lower extremities following the instructions in [53] were recorded using a 12-camera motion capture system (Vicon Motion Systems Ltd, Los Angeles, CA, USA) at 100 Hz. The GRF signals were collected at 1000 Hz synchronously with makers trajectories using in-ground force plates (AMTI, Watertown, MA, USA) mounted on an instrumented treadmill (Bertec Corp., Columbus, OH, USA) through the commercial real-time data capture software Nexus 2.9. Both GRF signals and markers trajectories were low-pass filtered with a fourth-order Butterworth filter in Visual 3D software (C-Motion, Rockville, MD, USA), and the cut-off frequencies were set as 25 Hz and 6 Hz, respectively. The markers trajectories and GRF signals from static poses and dynamic gait trials were used for joints kinematics and net joint moment calculation by using inverse kinematics and ID algorithms in Visual 3D. During all trials, four wired sEMG sensors (SX230, Biometrics Ltd, Newport, UK) were attached to the shank skin through a double-sided tape to non-invasively record the electrical signals from tibialis anterior (TA), LGS, MGS, and SOL muscles at 1000 Hz synchronously through Nexus 2.9. The locations for sEMG sensors were determined following the instructions in [54]. A US transducer (38-mm of length, 6.4 MHz center frequency, L7.5SC Prodigy Probe, S-Sharp, Taiwan) setup was applied in the walking experiments as described in [34]. The US radio frequency (RF) data were recorded at 1000 frames per second synchronously with markers trajectories, GRF, and sEMG signals, using a pulse sequence trigger signal generated from Nexus 2.9.

sEMG signals were band-pass filtered with the bandwidth between 20 Hz and 450 Hz, and then full-wave rectified and low-pass filtered with a cut-off frequency of 6 Hz. The resulting linear envelopes were normalized to the peak processed sEMG values obtained from all sets of dynamic gait trials with different speeds. The US radio frequency data were converted to B-mode images through beamforming and logarithmic compression, and then a commercial US imaging processing toolbox UltraTrack [55] was applied to determine the temporal sequences of the LGS and SOL muscles’ thickness based on the adaptive optical flow tracking algorithms. First, one region of interest (width by height: 336 \(\times\) 400 pixels) that encompassed the LGS and SOL muscles was selected as the area between the superficial and deep aponeuroses. Then, 10 vertical lines were manually defined with evenly distributed distances in the LGS’s region and SOL’s region on the first US imaging frame for each recorded waking trial, as shown in Fig. 2b. Key-frame correction [55] was applied to minimize the time-related drift of MT’s cyclical pattern over multiple gait cycles, where the key frames were selected to be at heel-strike and toe-off time points. After each correction, the new key frames’ vertical lines positions were determined by applying an affine transformation to the key-frame before it. Finally, the temporal sequence of the mean value of those 10 vertical lines’ lengths from LGS and SOL muscles were calculated for LGS’s MT and SOL’s MT, respectively. The MT signals were then low-pass filtered with a cut-off frequency of 30 Hz.

Before data collection, participants practiced dynamic walking steadily with all sensing devices attached to the lower extremities on the treadmill for at least 30 s for each speed. The duration of each walking trial was set as 2 min, and data from the middle 20 s were collected for analysis. In the current study, there were five different walking speeds, 0.50 m/s, 0.75 m/s, 1.00 m/s, 1.25 m/s, and 1.50 m/s, and the order was selected randomly for each participant. The participants were provided with at least 2 min for rest between two successive dynamic walking trials to avoid muscle fatigue.

Data post-processing

HNM model calibration

To determine the values for a set of parameters, including the shape factor \(A_{i}\), the tendon slack length \(l_{t_{i}}^{sk}\), the muscle contraction force at the walking task-specific maximum voluntary contraction \(F_{i}^{\max }\), and the muscle activation allocation gain \(\delta _{i}\), the HNM model calibration was performed with the initial parameters obtained from the literature [56]. The summarized calibration process of the sEMG-US imaging-driven HNM is presented in Fig. 2b. Indeed, the model calibration is a nonlinear optimization problem, where the above parameters set needs to be found for minimizing the designated objective function. During calibration, the parameters were subsequently adjusted within the predefined boundaries, which ensured the muscle-tendon unit always operated within the physiological range [39]. Except for the allocation gains \(\delta _{i}\) were set the range between 0 and 1, other parameters were set with the lower and upper bounds as 50% and 150% of the literature-referred values, while the initial guess of the set of parameters were set as the middle value of each parameter. Other HNM parameters, like muscle-tendon unit length \(l_{mt_{i}}\), optimal muscle fascicle length \(l_{mo_{i}}^{0}\), and optimal pennation angle \(\phi _{i}^{0}\) were assigned using OpenSim and the scaling method as references. The Matlab function ’lsqcurvefit’ with the Levenberg-Marquardt algorithm was applied to solve the nonlinear least-squares optimization problem until the following objective function was minimized during the calibration

$$\begin{aligned} E_{obj}=\frac{1}{N}\sum _{j=1}^{N}\left( \left( \sum _{i=1}^{2}M_{i}(j)\right) -M_{w}(j)\right) ^{2} \end{aligned}$$
(15)

where \(\sum _{i=1}^{2}M_{i}(j)\) represents the net PF moment estimation from both LGS and SOL muscles during stance phase at time instant j, \(M_{w}(j)\) represents the net PF moment measurement from ID at time instant j, and N represents the length of the data used for the calibration.

It should be noted that the HNM calibration is subjective, and for each participant, the calibration procedure was repeated for both single-speed modes and inter-speed mode. In the current study, single-speed modes are defined as the HNM calibration only with data collected from one single speed (five single-speed modes here), while the inter-speed mode is defined as the HNM calibration with all data collected from five different speeds (one inter-speed mode here). In HNM calibration, collected data, including ankle joint’s kinematics, kinetics, sEMG signals, and US imaging from the stance phase of 10 steady gait cycles, were used to minimize the objective function in (15). For the single-speed modes, the 10 steady gait cycles were all from one speed, while for the inter-speed mode, the 10 gait cycles were composed of two steady cycles from each of five speeds. Similarly, by manually setting the muscle activation allocation gain \(\delta _{i}\) to be 1 and 0, we could also conduct the calibration procedures of the sEMG-driven and US imaging-driven HNMs. After the HNM calibration with different modes, new data sets (five gait cycles not involved in the calibration procedures) from walking trials at different speeds were used as input of the calibrated HNMs to predict net PF moment. The prediction results were evaluated through the comparison to the measured net PF moment from ID.

Statistical analysis

The validation procedures comprised three tests to assess the sEMG-US imaging-driven, sEMG-driven, and US imaging-driven HNMs’ calibration and prediction abilities of the net ankle joint PF moment during the stance phase. In each type of HNM test, the root mean square error (RMSE), RMSE normalized to individual peak net PF moment (\(N-RMSE\)), RMSE normalized to individual body mass (\(BM-RMSE\)), and coefficient of determination (\(R^{2}\)) values between the calibrated/predicted net PF moment and ID-calculated net PF moment were calculated to evaluate the calibration/prediction performance. The rational of using these three tests is given here. From a mathematical perspective, the proposed sEMG-US imaging-driven HNM is a model-based regression problem. Given that the corresponding data were all collected continuously and synchronously, and the ground truth from the ID calculation and the prediction from different HNMs are all continuous, the most intuitive metric to evaluate the prediction accuracy is to use the RMSE between the ground truth and prediction. However, due to the weight, height, and walking pattern’s variations among subjects, the capability (peak value) of net plantarflexion moment generation during treadmill walking at the same speed is different, which would affect the prediction RMSE directly. Therefore, the direct RMSE values comparison among subjects is likely to introduce high deviations. An effective way to reduce the subjective variation of the prediction RMSE is to take the normalization to a subjective characteristic, like the body mass and peak net plantarflexion moment as we selected. The relationship between \(N-RMSE\) and \(BM-RMSE\) is that in the normalization calculation, they have same numerators but different denominators for the same subject. Although \(N-RMSE\) and \(BM-RMSE\) might provide potential redundancy information as can be seen from subsequent calibration and prediction summary results, they show different sensitivities to the walking speeds (\(BM-RMSE\) is more sensitive). Another useful metric to evaluate the prediction accuracy is the coefficient of determination, known as \(R^{2}\), which ranges from 0 to 1 and measures the proportion of variation in the data that is accounted for in the model. In the current study, \(R^{2}\) is used as a parallel evaluation metric as N-RMSE and BM-RMSE, and it is relatively not sensitive to the walking speeds. Although \(R^{2}\) and \(N-RMSE\) might provide potential redundancy information, it is not always true that lower \(N-RMSE\) corresponds to higher \(R^{2}\), so we keep both to provide potential supplementary information.

Shapiro-Wilk parametric hypothesis test was used to determine the normality of the corresponding \(N-RMSE\), \(BM-RMSE\), and \(R^{2}\) values of each calibration step under the single-speed modes and inter-speed mode, and each prediction step. According to the results of the Shapiro-Wilk test, two-way repeated-measure analysis of variance (ANOVA) or Friedman’s tests followed by a Tukey’s honestly significant difference tests (Tukey’s HSD) was applied to evaluate the effect of HNMs’ category and walking speed on different evaluation criteria, including \(N-RMSE\), \(BM-RMSE\), and \(R^{2}\) values, during both calibration and prediction procedures. The significant difference level was chosen as \(p<0.05\) for all statistical tests. Effect sizes were reported as \(\eta _{p}^{2}\) and Cohen’s d for main effects from ANOVA or Friedman’s tests and pairwise comparisons from Tukey’s HSD, respectively.

Results

Ankle joint kinematics and kinetics, and plantarflexors neuromuscular features

The treadmill walking speed’s effect on human ankle joint kinematics, kinematics, and lower extremities’ muscles activities has been extensively investigated and discussed in a recent study [57], but without results related to the architectural change of those muscles. In this study, the time sequences of both LGS and SOL’s muscle thickness during the recorded walking duration were extracted by using UltraTrack. Take the SOL muscle of Participant Sub08 as an example, Fig. 3 shows the MT tracking results throughout the 20 s walking duration at each walking speed. The red solid and blue dashed curves represent the MT change of the SOL muscle with and without key-frame correction, respectively. The results demonstrate that the key-frame correction could significantly reduce the time-related drift, which is due to the tracking error accumulation along with the walking duration (exampled US imaging video with the LGS and SOL’s MT tracking results at multiple walking speeds can be referred to Additional files 2, 3, 4, 5, 6). The reported results here are consistent with the studies mentioned in [55].

Fig. 3
figure 3

SOL muscle thickness tracking results by using UltraTrack at different walking speeds on Participant Sub08. Each subplot represents the tracking results during the recorded 20 s experiments under each walking speed, and the red solid and blue dashed lines represent the results with and without key-frame correction, respectively

The continuous results of the ankle joint kinematics and kinetics throughout the recorded 20 seconds at each speed were segmented as a percent of gait cycle from 0 to 100% according to the GRF measurements, which are shown in Fig. 4. The solid lines and light shadowed areas report the mean and one standard deviation (SD) values of the right ankle joint angular position, velocity, and net PF moment changes across all gait cycles at each walking speed on Sub01 (data from other participants can be found in the Additional file 1: Fig. S1–S5), where 0% and 100% represent the time instants when the heel-strike occurred in the current and consecutive gait cycles, respectively. By comparing the ankle joint moment calculated from ID, the averaged peak net PF moment across all walking gait cycles within each walking trial monotonically increased with the increase of walking speeds on each participant, indicating more power was needed for faster walking speed. The continuous results of neuromuscular features throughout the recorded walking duration were segmented as a percent of stance cycle from 0 to 100%, which are shown in Fig. 5. The solid lines and light shadowed areas report the mean and one SD values of the processed sEMG signals and US imaging-derived MT signals from the LGS and SOL muscles on the right leg during multiple stance phase cycles at each walking speed on Sub01 (data from other participants can be found in the Additional file 1: Fig. S11–S15), where 0% and 100% represent the time instants when the heel-strike and toe-off occurred in the same gait cycle, respectively. The upper two rows present the LGS and SOL MT changes during the stance phase, while the lower two rows present the muscle’s low-pass filtered sEMG changes. The results show small variations of MT and sEMG for both LGS and SOL muscles at heel-strike across all walking speeds, but relatively high variations of sEMG for both LGS and SOL muscles at toe-off, where the sEMG signals increase with increased walking speed. For all five walking speeds, we observe that MT of LGS or SOL muscle is almost identical at heel-strike and toe-off, however, the processed sEMG signals of LGS or SOL muscle show higher value at toe-off than at heel-strike. By comparing the shadowed areas of the same feature across speeds, we observed that features’ deviations are higher at a slower speed, like 0.50 m/s, than those at higher speeds. This implies that keeping a steady muscle contraction pattern at a slower speed is more difficult than at a higher speed.

Fig. 4
figure 4

Results of joint kinematics and kinetics during the treadmill walking experiments at different speeds on Participant Sub01. The mean and one SD values on the right ankle joint during recorded gait cycles (between the heel-strike instants in the current and the consecutive gait cycles) are represented by the solid lines and the light shadowed areas, respectively

Fig. 5
figure 5

Results of processed sEMG signals and US imaging-derived MT signals from both LGS and SOL muscles during the treadmill walking experiments at different speeds on Participant Sub01. The mean and one SD values of the plantarflexor muscles on the right shank during recorded stance phase cycles (from the heel-strike instant to the toe-off instant in the same gait cycle) are represented by the solid lines and the light shadowed areas, respectively

Correlation analysis results

Overall, in Fig. 5, both temporal MT and sEMG signals from LGS and SOL muscles show a fairly strong correlation with the net ankle joint PF moment calculated from ID. For each participant, it was assumed the step-to-step variation is negligible, so we segmented all periods of the ankle joint net plantarflexion moment during the walking stance phase at each walking speed, and then calculated a correlation coefficient value between the net plantarflexion moment time sequence and each neuromuscular feature time sequence during the same stance cycle. Therefore, for a specific neuromuscular feature and walking speed, if the stance cycles are \(n_{1}, n_{2}, ..., n_{10}\) for all ten participants, the number of correlation coefficients would be \(\sum _{i=1}^{10}n_{i}\). To address the person-to-person variation and to get each entry in Table 1, we applied the Fisher transformation to yield approximately normally distributed correlation coefficients. After the Fisher transformation, we calculated the mean and one SD values of the correlation coefficients for each walking speed and the neuromuscular feature, as listed in Table 1. The results show that a positive relationship exists between the processed LGS or SOL’s sEMG signal and net PF moment, and between the LGS or SOL’s MT and net PF moment. All the mean values of correlation coefficients are higher than 0.8 except for the SOL’s sEMG signal at 1.50 m/s walking speed. From the results of a two-way repeated-measure ANOVA, we did not observe a significant difference in the correlation coefficient values among the four neuromuscular features (\(p=0.201\)) or among the walking speed (\(p=0.112\)), which indicates the MT and sEMG signals might have comparable capabilities to predict net PF moment across speeds.

Table 1 Correlation coefficients (mean ± SD) between each neuromuscular feature and the net PF moment across all stance phase cycles and ten participants at each walking speed

Results of HNM calibration and net PF moment prediction

Fig. 6
figure 6

Muscle activation levels of both LGS and SOL during the walking stance phase by setting the allocation gain as 0, 1, and the optimal value, respectively. The blue solid, red dashed, and orange centered curves represent the muscle activation levels by using only sEMG signals, only US imaging-derived MT signals, and the fusion between sEMG and MT signals with the optimal allocation gain in both model calibration (left column) and prediction (right column) procedures. Data shown in the figure are from Participant Sub03 when walking on the treadmill at 0.75 m/s

In the HNM calibration procedures shown in Fig. 2b, one essential input component is the muscle activation levels of both LGS and SOL muscles, either only based on sEMG signals, US imaging-derived MT signals, or data fusion between them. The representative demonstration of three types of muscle activation levels of both LGS and SOL muscles during the walking stance phase at 0.75 m/s on Participant Sub03 are shown in Fig. 6. The left and right column subplots represent the time sequence points of muscle activation levels that are used in the model calibration and prediction procedures. It is observed that the data fusion could balance the muscle activation levels from both sEMG signals and MT signals and further compensate for the activation level drift from only MT signals, which is potentially beneficial for the accuracy improvement of the net plantarflexion moment prediction. Fig. 7 presents the representative results of HNMs calibrations on Sub05 with data collected from each individual walking speed, where each solid centered line and the shadowed area represent the mean and one SD values of the measured and calibrated net ankle joint PF moment, and each curve is composed of 10 stance cycles at each walking speed. For this participant, the mean calibration RMSE values are 8.70 N\(\cdot\)m, 9.10 N\(\cdot\)m, and 7.90 N\(\cdot\)m by using sEMG-, US imaging-, and sEMG-US imaging-driven HNMs under walking speed of 0.50 m/s, respectively. Results from other walking speeds also showed that the sEMG-US imaging-driven HNM could effectively reduce the calibration RMSE under either single-speed modes or inter-speed mode. To quantitatively evaluate the calibration performance by using different HNM categories and under different speed modes, the mean RMSE and \(R^{2}\) values between the benchmark moment values from ID and calibrated moment values on each participant are summarized in Table 2. Based on the observation from Fig. 7, it appears that the relative error may be inflated during specific areas of the stance cycles. The current outcomes (i.e., \(N-RMSE\), \(BM-RMSE\), and \(R^2\) values) may not be directly sensitive to those specific areas, because the three evaluation metrics report an averaged performance across the entire stance cycle. This limitation of the current outcome metrics can be addressed in future work by identifying systematic changes in the prediction accuracy throughout the stance phase cycle rather than across the entire cycle.

Fig. 7
figure 7

Calibration results by using three HNM categories on Participant Sub05 under single-speed modes. The calibration performance by using sEMG-, US imaging-, and sEMG-US imaging-driven HNMs under 5 single-speed modes on Sub05. The red lines and light shadowed areas represent the mean and SD of the net PF moment from the ID calculation, the blue, green, black lines and light shadowed areas represent the mean and SD of the net PF moment from calibrations by using sEMG-driven HNM, US imaging-driven HNM, and sEMG-US imaging-driven HNM, respectively. a–e represent the results at each walking speed from 0.50 to 1.50 m/s

Table 2 Calibration results with three HNM categories in both single-speed modes and inter-speed mode on each participant
Fig. 8
figure 8

Calibration performance across all participants by using three HNM categories. The mean and SD of \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\) values across all participants under both single-speed modes and inter-speed mode are displayed in upper, middle, and lower plots, respectively, during calibration procedures. Blue, red, and yellow bars represent the results by using sEMG-, US imaging-, and sEMG-US imaging-driven HNMs, respectively. Letters and asterisks represent the statistically significant difference from the two-way repeated-measure ANOVA followed by a Tukey’s HSD, where letters stand for the evaluation of the walking speed factor (significant difference only exists between each two bar groups out of the left six groups without any overlap of letters) and asterisks stand for the evaluation of the applied HNM type factor (significant difference exists between each two bars out of the most right group when the asterisk is there)

Overall, the calibration RMSE/\(R^{2}\) values by using the sEMG-US imaging-driven HNM are lower/higher than those by using the sEMG-driven and US imaging-driven HNMs, respectively, for both single-speed modes and inter-speed mode on all participants. Additionally, 97.5 % of \(R^{2}\) values are higher than 0.80, which indicates a strong linear correlation between the calibrated and ID-calculated net PF moment values. After calculating the \(N-RMSE\) and \(BM-RMSE\) values, the summarized calibration performance across all participants under both single-speed modes and inter-speed mode is presented in Fig. 8. Shapiro-Wilk tests show that all \(N-RMSE\), \(BM-RMSE\), and \(R^{2}\) values under the single-speed modes and inter-speed mode in both calibration and prediction are normally distributed. In calibration, results across participants from the two-way ANOVA show \(BM-RMSE\) values systematically changed in response to changes in speed (main effect, \(p<0.001\), \(n_{p}^{2}=0.268\)) and applied HNMs (main effect, \(p<0.001\), \(n_{p}^{2}=0.162\)). \(N-RMSE\) values systematically changed in response to changes in speed (main effect, \(p<0.001\), \(n_{p}^{2}=0.089\)) and applied HNMs (main effect, \(p<0.001\), \(n_{p}^{2}=0.193\)). \(R^{2}\) values were significantly affected by speed (main effect, \(p<0.001\), \(n_{p}^{2}=0.153\)) and applied HNMs (main effect, \(p<0.001\), \(n_{p}^{2}=0.186\)). However, there is no evidence of an interaction effect of the speed and HNM on \(BM-RMSE\) values (\(p=0.679\)), \(N-RMSE\) values (\(p=0.948\)), and \(R^{2}\) values (\(p=0.279\)), respectively. Additionally, the results from the post hoc Tukey’s HSD are also marked in Fig. 8, where letters and asterisks represent the statistically significant difference.

To verify the net PF moment prediction performance of the proposed sEMG-US imaging-driven HNM and compare the new HNM to both sEMG- and US imaging-driven HNMs, the three calibrated HNM categories, shown in Fig. 8, were applied to walking scenarios of five different speeds. Similarly, \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\) values between the predicted and ID-calculated net PF moments were evaluated and compared. Due to the space limitation, the detailed RMSE and \(R^{2}\) values at each walking speed on individual participant in prediction can be found in Additional file 1: Tables S2 to Table S4. The summarized results in Fig. 9 show the mean and one SD of \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\) values in prediction at each of five speeds by applying three HNM categories. The subplots from the first column to the fifth column represent the applied calibrated HNMs were from five single-speed modes, while the subplots from the sixth column represent the inter-speed mode. In each subplot, the \(x-\)axis labels represent the prediction scenarios at five different speeds. Table 3 summaries the mean values of \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\) in the prediction by using the proposed sEMG-US imaging-driven HNM, where the bold numbers represent the cases that the calibration speed mode is the same as the prediction speed scenario, which always results in the lowest RMSE and highest \(R^{2}\) values. From Fig. 9 and Table 3, we observed that \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\) values are optimal when the applied calibrated single-speed mode is close to the current prediction speed scenario, regardless of HNM categories. Furthermore, the calibrated inter-speed mode can effectively constrain the \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\) values within a small variation range among all prediction speed scenarios, regardless of HNM categories.

The results of two-way repeated-measure ANOVA show that in prediction \(BM-RMSE\) values are significantly affected by the calibration speed mode (main effect, \(p<0.001\), \(\eta _{p}^{2}=0.104\)) and applied HNM category (main effect, \(p<0.001\), \(\eta _{p}^{2}=0.141\)), \(N-RMSE\) values are significantly affected by the calibration speed mode (main effect, \(p<0.001\), \(\eta _{p}^{2}=0.101\)) and applied HNM category (main effect, \(p<0.001\), \(\eta _{p}^{2}=0.164\)), and also \(R^{2}\) values are significantly affected by the calibration speed mode (main effect, \(p<0.001\), \(\eta _{p}^{2}=0.103\)) and applied HNM category (main effect, \(p<0.001\), \(\eta _{p}^{2}=0.164\)). The pairwise comparisons show that across all prediction scenarios with every speed mode calibration, the proposed sEMG-US imaging-derived HNM could significantly reduce the \(BM-RMSE\) values by 14.49% (\(p=0.029\), \(d=-0.348\)) and 36.94% (\(p<0.001\), \(d=-0.930\)), significantly reduce the \(N-RMSE\) values by 14.58% (\(p=0.012\), \(d=-0.389\)) and 36.79% (\(p<0.001\), \(d=-1.021\)), and increase the \(R^{2}\) values by 9.71% (\(p=0.023\), \(d=0.713\)) and 17.86% (\(p<0.001\), \(d=1.010\)) compared to sEMG-driven and US imaging-driven HNMs, respectively. Furthermore, the pairwise comparisons also show that across all prediction scenarios with every HNM category, the inter-speed mode calibration could reduce the \(BM-RMSE\) values by 38.74% (\(p<0.001\), \(d=-0.922\)), 18.86% (\(p=0.180\), \(d=-0.590\)), 8.87% (\(p=0.826\), \(d=-0.237\)), 9.84% (\(p=0.785\), \(d=-0.299\)), and 26.57% (\(p=0.004\), \(d=-0.836\)), reduce the \(N-RMSE\) values by 35.73% (\(p<0.001\), \(d=-1.002\)), 17.20% (\(p=0.201\), \(d=-0.403\)), 9.11% (\(p=0.881\), \(d=-0.254\)), 11.44% (\(p=0.715\), \(d=-0.326\)), and 28.08% (\(p<0.001\), \(d=-0.817\)), and increase the \(R^{2}\) values by 16.76% (\(p<0.001\), \(d=0.810\)), 5.33% (\(p=0.307\), \(d=0.342\)), 0.40% (\(p=0.896\), \(d=0.044\)), 0.31% (\(p=0.921\), \(d=0.033\)), and 6.51% (\(p=0.127\), \(d=0.487\)) compared to single-speed mode calibration at 0.50 m/s, 0.75 m/s, 1.00 m/s, 1.25 m/s, and 1.50 m/s, respectively.

Fig. 9
figure 9

Prediction performance across all participants by using three HNM categories. The mean and SD of \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\) values across all participants under both single-speed modes and inter-speed model are displayed in upper, middle, and lower row plots, respectively. Blue, red, and yellow bars represent the prediction results by using sEMG-, US imaging-, and sEMG-US imaging-driven HNMs, respectively. Column subplots from left to right represent that these applied HNMs were calibrated under single-speed modes (0.50 m/s, 0.75 m/s, 1.00 m/s, 1.25 m/s, and 1.50 m/s) and the inter-speed mode

Table 3 Metrics of \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\) mean values of the net PF moment prediction across all participants by using the proposed sEMG-US imaging-driven HNM

Discussions

For the first time, this study investigated net ankle joint PF moment continuous prediction during the walking stance phase based on an HNM approach that combines both sEMG-induced and US imaging-induced muscle activation components. We defined this approach as sEMG-US imaging-driven HNM for the net PF moment prediction. The modified HNM was validated with quasi-periodic data, including sEMG signals and US imaging-derived MT from plantarflexors, ankle joint kinematics, and kinetics measurements from multiple treadmill walking stance phases gait cycles at five different speeds. We focused on two plantarflexor muscles, e.g., LGS and SOL, when characterizing the PF function during walking. The basic objective was to establish a modified HNM for either LGS or SOL muscles, including (1) a weighted muscle activation model component, (2) a muscle-tendon unit geometry model component, and (3) a muscle dynamic contraction model component. This study also investigated the effects of HNM categories (sEMG-, US imaging-, and sEMG-US imaging-driven) and calibration speed modes (single-speed modes and inter-speed mode) on the net PF moment prediction performance during the walking stance phase. The ID-calculated net PF moment worked as the benchmark in both calibration and prediction. The results support our hypothesis that compared to sEMG-driven and US imaging-driven HNMs, the sEMG-US imaging-driven HNM accurately predicts the net PF moment. Further, the HNMs, when calibrated with the inter-speed mode, robustly predict the net PF moment.

The evaluation of the proposed HNMs’ calibration and prediction performance was based on three criteria: \(BM-RMSE\), \(N-RMSE\), and \(R^{2}\). In general, as reported in [58], the results were considered excellent if the \(N-RMSE\) values were smaller than 15%. For the calibration, the results in Fig. 7 presented a good calibration performance by using each of three HNMs on Sub05. Furthermore, from results across participants shown in Fig. 8, we observed that the mean \(N-RMSE\) values were all less than 10.92% regardless of the calibration speed modes and applied HNM categories. Primarily they were less than 8.30% irrespective of the calibration speed modes when applying the proposed sEMG-US imaging-driven HNM. For the prediction, as shown in Fig. 9 and Table 3, by applying the inter-speed mode calibrated sEMG-US imaging-driven HNM, the \(N-RMSE\) mean values were all less than 11.40% throughout those five walking speeds scenarios, which validated the prediction results were all excellent. Therefore, both calibration and prediction results present consistent findings that the \(N-RMSE\) values can be constrained within a satisfactorily small region by using the proposed sEMG-US imaging-driven HNM with the inter-speed mode calibration.

The superior net PF torque prediction of the sEMG-US imaging-driven HNM can be mainly due to the reason that the fused signals capture complementary mechanical and neural aspects of the muscle contraction [25]. information [25]. Specifically, sEMG signals measure electric potentials generated by muscle motor units. The amplitude and density of sEMG signals linearly correlate with the number of firing neurons, which offers a physical measurement of the microphysiological response [59]. US imaging directly visualizes muscle change at the macrophysiological performance [60], when the same group of neural motor units is activated. Thus, sEMG features and US imaging features provide information from electrical and mechanical aspects, respectively, in response to the same physiological stimulus.

As a preliminary study of the sEMG-US imaging-driven HNM for walking on a treadmill across multiple speeds, the results are promising and can help overcome the challenges of joint motion intent detection in volitional control of assistive devices. However, there are still some limitations in the current study. First of all, only a few US imaging-derived parameters, MT from both LGS and SOL, measured muscle activation levels indirectly (a macrophysiological perspective). Other US-derived parameters, like FL and PA [21, 25, 27, 35], of the MTU geometry model may further enrich the neuromuscular model. A potential difficulty here is the inability to capture a larger region of interest that contains entire muscle fascicles, mainly due to the small dimension of our US transducer. Another difficulty is to continuously track the muscle fascicles of the plantarflexors within the visualized region of interest during dynamic tasks, like walking. If always visualized, more geometry parameters from US imaging signals will likely further improve the HNM accuracy. Second, the US imaging-derived muscle activation was defined as a normalization function in (4). We also assumed a linear monotonic relationship between the MT of LGS or SOL and net PF moment change during the stance phase in Table 1. However, in [29, 61], it was observed that despite an increase in MT of the LGS and SOL muscles during the stance phase, their muscle activations decreased in the last 10 to 20% of the stance phase. Also, in fast walking, MT of SOL was observed to decrease in the first 10% of stance, while the muscle activation increased. The above inconsistency may result from the difference in experimental setup, sEMG sensor or US transducer placement on the targeted muscles, MT tracking approach from US imaging signals, etc. In future work, a further comparison will be conducted to evaluate the muscle geometry parameters change and sEMG signals change during the same walking tasks. As detailed in the pre-processing data section, a key-frame correction [55] was applied to reduce the gait cycle-related drift of MT tracking results significantly. However, a subjective error in the manual operation is inevitable. Therefore, it was not practical to altogether remove the gait cycle-related drift of US imaging-derived MT. This work has extended the preliminary results of voluntary isometric ankle PF studies [34, 35] to dynamic walking studies on non-disabled participants. Further experiments on patients with impaired plantarflexors are needed to verify the practicability of the proposed sEMG-US imaging-driven HNM for net ankle joint PF moment prediction. We also observed that the calibrated HNMs with the inter-speed mode increased the net PF moment prediction robustness across different speed scenarios. The development of inter-speed HNM modeling is crucial because it will potentially allow a consistently accurate prediction performance of the net PF moment across multiple walking speeds. Furthermore, it will enable predicting and analyzing the dynamics of a more extensive set of walking speeds than was possible with a single-speed mode, resulting in a deeper understanding of the neuromuscular dynamics during the human walking stance phase. Finally, it may open up to the development of robust neuromuscular human-machine interfaces for the simultaneous and proportional control of wearable assistive devices such as powered orthoses and prostheses under different walking speeds.

In this preliminary study, although the combination of US imaging and sEMG data in the model improved the prediction of the ankle moment, it was unknown if the prediction was better only because in this combination the solution space was higher by adding parameters in the optimization problem. To address this ambiguity, we have conducted additional model calibration and prediction procedures for comparison, where the time sequence of the US imaging-derived muscle activation was replaced by a randomly generated signal between 0 and 1. In this way, we could reproduce the same solution space just like the proposed sEMG-US imaging-driven HNM, so we would be able to investigate if the prediction benefits still exist or not. We took a representative participant as an example, and present the calibration performance and prediction performance by using five different HNMs, as shown in the Additional file 1: Fig. S18. It was observed that when the random signal (RS)-derived muscle activation levels of LGS and SOL muscles were introduced, although the solution space was increased, just like the proposed approach to fuse sEMG and US imaging signals, the calibration and prediction performance were both deteriorated when compared to the case where only sEMG-derived muscle activation levels were used. Therefore, based on the additional analyses, we believe that the better net plantarflexion moment prediction by using the proposed sEMG-US imaging-driven HNM was not only because in this combination the solution space was higher by adding parameters in the optimization problem, but more importantly, the US imaging-derived muscle thickness feature introduced complementary mechanical deformation characteristics during the muscle contraction, which dominantly resulted in the better prediction performance.

Limitation and future work

The key point to implement the proposed sEMG-US imaging-driven HNM is to find the allocation gain \(\delta _{i}\) for both LGS and SOL muscles through solving the nonlinear least square optimization problem. However, one possible limitation of this optimization process is that the allocation gains are speed-dependent and subject-dependent, which means that the optimal values of those two allocation gains vary along with the model calibrations when using data collected from different walking speeds on the same participant or from different participants at the same walking speed. So it is relatively challenging to generate a constant allocation gain for either LGS or SOL muscles across multiple speeds or general participants. The US imaging-derived muscle thickness signal is potentially sensitive to the position and orientation of the US transducer relative to the shank segment, which may affect the prediction accuracy of the proposed approach. In the current study, we did not track the position and orientation of the US transducer relative to the shank segment due to the assumption that the relative motion of the US transducer is negligible because a customized US transducer holder with an arc structure and elastic velcro straps were used to stabilize the position and orientation of the US transducer during the dynamic walking tasks. However, one possible limitation of this study is that the assumption may not hold for a faster walking speed. Therefore, a further study would be to investigate the effect of the relative position and orientation change between the US transducer and the shank segment to the prediction accuracy of the proposed sEMG-US imaging-driven HNM approach.

Another limitation of the current study is that we only considered the agonist LGS and SOL muscles to predict the ankle joint net PF moment, and the effects of the MGS and TA muscles were not considered. Although in the experiments we collected the sEMG signals from LGS, MGS, SOL, and TA muscles, ultrasound imaging data was collected from the LGS and SOL muscles only. Therefore, the fusion results were only based on these two muscles. The main reasons that we did not consider the effects from the MGS and TA muscles are (1) existing evidence suggests that sEMG signals from LGS and MGS are comparable [62], and similar results in the current study also can be seen in the Additional file 1: Fig. S19; (2) during the walking stance phase, plantarflexor muscles (LGS and SOL) are the dominant muscles that generate the net ankle joint moment. Results of the processed sEMG signals showed that the activation of the TA muscle are relatively small during the walking stance phase compared to the plantarflexor muscles, as can be observed in the Additional file 1: Fig. S19. However, even though the TA muscle activation is low during the walking stance phase, the effect from the passive contraction from this antagonist muscle needs to be further investigated in future work.

Conclusion

This paper investigated the feasibility of an sEMG-US imaging-driven HNM to predict net ankle joint PF moment during the stance phase across multiple walking speeds. The results showed that on average, the net PF moment prediction RMSE, normalized to peak net PF moment/body mass, was significantly reduced when using the sEMG-US imaging-driven HNM, compared to sEMG-driven and US imaging-driven HNMs. The results also showed that the calibrated HNMs with the inter-speed mode were more robust for net PF moment prediction at different speed scenarios than HNMs with single-speed modes. The improved net ankle joint PF moment prediction during the dynamic walking tasks at different speeds could potentially lead to improved volitional control of assistive devices with more advanced and intelligent algorithms.

Availability of data and materials

The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.

Abbreviations

GRF:

Ground reaction force

ID:

inverse dynamics

HNM:

Hill-type neuromuscular model

EMG:

Electromyography

sEMG:

Surface electromyography

PA:

Pennation angle

FL:

Fascicle length

MT:

Muscle thickness

MGS:

Medial gastrocnemius

LGS:

Lateral gastrocnemius

SOL:

Soleus

PF:

Plantarflexion

EMD:

Electromechanical delay

MTU:

Muscle-tendon unit

IRB:

Institutional Review Board

B-mode:

Brightness mode

ANOVA:

Analysis of variance

HSD:

Honestly significant difference test

RMSE:

Root mean square error

N-RMSE:

Root mean square error normalized to peak value

BM-RMSE:

Root mean square error normalized to body mass

\({\eta _p}^2\) :

Effect size

SD:

Standard deviation

R :

Correlation coefficient

\(R^2\) :

Coefficient of determination

d :

Cohen’s d

US:

Ultrasound

References

  1. Tzu-wei PH, Shorter KA, Adamczyk PG, Kuo AD. Mechanical and energetic consequences of reduced ankle plantar-flexion in human walking. J Exp Biol. 2015;218(22):3541–50.

    Google Scholar 

  2. Iezzoni LI, McCarthy EP, Davis RB, Siebens H. Mobility difficulties are not only a problem of old age. J Gen Intern Med. 2001;16(4):235–43.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  3. Takahashi KZ, Lewek MD, Sawicki GS. A neuromechanics-based powered ankle exoskeleton to assist walking post-stroke: a feasibility study. J Neuroeng Rehabil. 2015;12(1):23.

    PubMed  PubMed Central  Article  Google Scholar 

  4. Tamburella F, Tagliamonte N, Pisotta I, Masciullo M, Arquilla M, Van Asseldonk E, van Der Kooij H, Wu A, Dzeladini F, Ijspeert A. Neuromuscular controller embedded in a powered ankle exoskeleton: Effects on gait, clinical features and subjective perspective of incomplete spinal cord injured subjects. IEEE Trans Neural Syst Rehabil Eng. 2020;28(5):1157–67.

    CAS  PubMed  Article  Google Scholar 

  5. Young AJ, Ferris DP. State of the art and future directions for lower limb robotic exoskeletons. IEEE Trans Neural Sys Rehabil Eng. 2017;25(2):171–82.

    Article  Google Scholar 

  6. Quinlivan B, Lee S, Malcolm P, Rossi D, Grimmer M, Siviy C, Karavas N, Wagner D, Asbeck A, Galiana I. Assistance magnitude versus metabolic cost reductions for a tethered multiarticular soft exosuit. Sci Robot. 2017;2(2):4416.

    Article  Google Scholar 

  7. Asbeck AT, De Rossi SM, Holt KG, Walsh CJ. A biologically inspired soft exosuit for walking assistance. Int J Robot Res. 2015;34(6):744–62.

    Article  Google Scholar 

  8. Kesar TM, Perumal R, Reisman DS, Jancosko A, Rudolph KS, Higginson JS, Binder-Macleod SA. Functional electrical stimulation of ankle plantarflexor and dorsiflexor muscles: effects on poststroke gait. Stroke. 2009;40(12):3821–7.

    PubMed  PubMed Central  Article  Google Scholar 

  9. Embrey DG, Holtz SL, Alon G, Brandsma BA, McCoy SW. Functional electrical stimulation to dorsiflexors and plantar flexors during gait to improve walking in adults with chronic hemiplegia. Arch Phys Med Rehabil. 2010;91(5):687–96.

    PubMed  Article  Google Scholar 

  10. Sharma N, Stegath K, Gregory CM, Dixon WE. Nonlinear neuromuscular electrical stimulation tracking control of a human limb. IEEE Trans Neural Syst Rehabil Eng. 2009;17(6):576–84.

    PubMed  Article  Google Scholar 

  11. Sharma N, Gregory C, Dixon WE. Predictor-based compensation for electromechanical delay during neuromuscular electrical stimulation. IEEE Trans Neural Syst Rehabil Eng. 2011;19(6):601–11.

    PubMed  Article  Google Scholar 

  12. Alibeji N, Kirsch N, Farrokhi S, Sharma N. Further results on predictor-based control of neuromuscular electrical stimulation. IEEE Trans Neural Syst Rehabil Eng. 2015;23(6):1095–105.

  13. Hussain S, Jamwal PK, Ghayesh MH, Xie SQ. Assist-as-needed control of an intrinsically compliant robotic gait training orthosis. IEEE Trans Ind Electron. 2016;64(2):1675–85.

    Article  Google Scholar 

  14. Selfslagh A, Shokur S, Campos DS, Donati AR, Almeida S, Yamauti SY, Coelho DB, Bouri M, Nicolelis MA. Non-invasive, brain-controlled functional electrical stimulation for locomotion rehabilitation in individuals with paraplegia. Sci Rep. 2019;9(1):1–17.

    CAS  Article  Google Scholar 

  15. Wang S, Wang L, Meijneke C, van Asseldonk E, Hoellinger T, Cheron G, Ivanenko Y, La Scaleia V, Sylos-Labini F, Molinari M, Tamburella F, Pisotta I, Thorsteinsson F, Ilzkovitz M, Gancet J, Nevatia Y, Hauffe R, Zanow F, van der Kooij H. Design and control of the MINDWALKER exoskeleton. IEEE Trans Neural Syst Rehabil Eng. 2015;23(2):277–86.

    PubMed  Article  Google Scholar 

  16. Zanotto D, Akiyama Y, Stegall P, Agrawal SK. Knee joint misalignment in exoskeletons for the lower extremities: effects on user’s gait. IEEE Trans Robot. 2015;31(4):978–87.

    Article  Google Scholar 

  17. Sartori M, Farina D, Lloyd DG. Hybrid neuromusculoskeletal modeling to best track joint moments using a balance between muscle excitations derived from electromyograms and optimization. J Biomech. 2014;47(15):3613–21.

    PubMed  Article  Google Scholar 

  18. Zhang Q, Sheng Z, Moore-Clingenpeel F, Kim K, Sharma N. Ankle dorsiflexion strength monitoring by combining sonomyography and electromyography. In: 2019 IEEE 16th International Conference on Rehabilitation Robotics (ICORR), 2019;pp. 240–245. IEEE

  19. Han J, Ding Q, Xiong A, Zhao X. A state-space emg model for the estimation of continuous joint movements. IEEE Trans Ind Electron. 2015;62(7):4267–75.

    Article  Google Scholar 

  20. Batzianoulis I, Krausz NE, Simon AM, Hargrove L, Billard A. Decoding the grasping intention from electromyography during reaching motions. J NeuroEng Rehabil. 2018;15(1):57.

    PubMed  PubMed Central  Article  Google Scholar 

  21. Zhang Q, Iyer A, Sun Z, Kim K, Sharma N. A dual-modal approach using electromyography and sonomyography improves prediction of dynamic ankle dorsiflexion motion. IEEE Trans Neural Syst Rehabil Eng. 2021;29:1944–54.

    PubMed  Article  Google Scholar 

  22. Hargrove LJ, Englehart K, Hudgins B. A comparison of surface and intramuscular myoelectric signal classification. IEEE Trans Biomed Eng. 2007;54(5):847–53.

    PubMed  Article  Google Scholar 

  23. Crouch DL, Pan L, Filer W, Stallings JW, Huang H. Comparing surface and intramuscular electromyography for simultaneous and proportional control based on a musculoskeletal model: a pilot study. IEEE Trans Neural Syst Rehabil Eng. 2018;26(9):1735–44.

    PubMed  Article  Google Scholar 

  24. Strasser EM, Draskovits T, Praschak M, Quittan M, Graf A. Association between ultrasound measurements of muscle thickness, pennation angle, echogenicity and skeletal muscle strength in the elderly. Age. 2013;35(6):2377–88.

    PubMed  PubMed Central  Article  Google Scholar 

  25. Zhang Q, Kim K, Sharma N. Prediction of ankle dorsiflexion moment by combined ultrasound sonography and electromyography. IEEE Trans Neural Syst Rehabil Eng. 2020;28(1):318–27.

    PubMed  Article  Google Scholar 

  26. Arampatzis A, Karamanidis K, Stafilidis S, Morey-Klapsing G, DeMonte G, Brüggemann G-P. Effect of different ankle-and knee-joint positions on gastrocnemius medialis fascicle length and emg activity during isometric plantar flexion. J Biomech. 2006;39(10):1891–902.

    PubMed  Article  Google Scholar 

  27. Zhang Q, Iyer A, Kim K, Sharma N. Evaluation of non-invasive ankle joint effort prediction methods for use in neurorehabilitation using electromyography and ultrasound imaging. IEEE Trans Biomed Eng. 2021;68(3):1044–55.

    PubMed  Article  Google Scholar 

  28. Damiano DL, Prosser LA, Curatalo LA, Alter KE. Muscle plasticity and ankle control after repetitive use of a functional electrical stimulation device for foot drop in cerebral palsy. Neurorehabil Neural Repair. 2013;27(3):200–7.

    PubMed  Article  Google Scholar 

  29. Hodson-Tole E, Lai A. Ultrasound-derived changes in thickness of human ankle plantar flexor muscles during walking and running are not homogeneous along the muscle mid-belly region. Sci Rep. 2019;9(1):1–11.

    CAS  Article  Google Scholar 

  30. Guo J-Y, Zheng Y-P, Xie H-B, Chen X. Continuous monitoring of electromyography (EMG), mechanomyography (MMG), sonomyography (SMG) and torque output during ramp and step isometric contractions. Med Eng Phys. 2010;32(9):1032–42.

    PubMed  Article  Google Scholar 

  31. Shi J, Zheng YP, Huang QH, Chen X. Continuous monitoring of sonomyography, electromyography and torque generated by normal upper arm muscles during isometric contraction: sonomyography assessment for arm muscles. IEEE Trans Biomed Eng. 2008;55(3):1191–8.

    PubMed  Article  Google Scholar 

  32. Dick TJM, Biewener AA, Wakeling JM. Comparison of human gastrocnemius forces predicted by Hill-type muscle models and estimated from ultrasound images. J Exp Biol. 2017;220(9):1643–53.

    PubMed  PubMed Central  Google Scholar 

  33. Jahanandish MH, Fey NP, Hoyt K. Lower-limb motion estimation using ultrasound imaging: a framework for assistive device control. IEEE J Biomed Health Inform. 2019;23(6):2505–2514

  34. Zhang Q, Iyer A, Kim K, Sharma N. Volitional contractility assessment of plantar flexors by using non-invasive neuromuscular measurements. In: 2020 8th IEEE RAS/EMBS International Conference for Biomedical Robotics and Biomechatronics (BioRob), pp. 515–20. IEEE

  35. Zhang Q, Clark WH, Franz JR, Sharma N. Personalized fusion of ultrasound and electromyography-derived neuromuscular features increases prediction accuracy of ankle moment during plantarflexion. Biomed Signal Process Control. 2022;71: 103100.

    Article  Google Scholar 

  36. Sikdar S, Rangwala H, Eastlake EB, Hunt IA, Nelson AJ, Devanathan J, Shin A, Pancrazio JJ. Novel method for predicting dexterous individual finger movements by imaging muscle activity using a wearable ultrasonic system. IEEE Trans Neural Syst Rehabil Eng. 2013;22(1):69–76.

    PubMed  Article  Google Scholar 

  37. Dhawan AS, Mukherjee B, Patwardhan S, Akhlaghi N, Diao G, Levay G, Holley R, Joiner WM, Harris-Love M, Sikdar S. Proprioceptive sonomyographic control: a novel method for intuitive and proportional control of multiple degrees-of-freedom for individuals with upper extremity limb loss. Sci Rep. 2019;9(1):1–15.

    Article  CAS  Google Scholar 

  38. Sheng Z, Sharma N, Kim K. Quantitative assessment of changes in muscle contractility due to fatigue during nmes: an ultrasound imaging approach. IEEE Trans Biomed Eng. 2019;67(3):832–41.

    PubMed  Article  Google Scholar 

  39. Lloyd DG, Besier TF. An EMG-driven musculoskeletal model to estimate muscle forces and knee joint moments in vivo. J Biomech. 2003;36(6):765–76.

    PubMed  Article  Google Scholar 

  40. Hodges PW, Pengel LHM, Herbert RD, Gandevia SC. Measurement of muscle contraction with ultrasound imaging. Muscle Nerve. 2003;27(6):682–92.

    CAS  PubMed  Article  Google Scholar 

  41. Shi J, Zheng Y-P, Chen X, Huang Q-H. Assessment of muscle fatigue using sonomyography: muscle thickness change detected from ultrasound images. Med Eng Phys. 2007;29(4):472–9.

    CAS  PubMed  Article  Google Scholar 

  42. McMeeken J, Beith I, Newham D, Milligan P, Critchley D. The relationship between emg and change in thickness of transversus abdominis. Clin Biomech. 2004;19(4):337–42.

    CAS  Article  Google Scholar 

  43. Kiesel KB, Uhl TL, Underwood FB, Rodd DW, Nitz AJ. Measurement of lumbar multifidus muscle contraction with rehabilitative ultrasound imaging. Man Ther. 2007;12(2):161–6.

    PubMed  Article  Google Scholar 

  44. Zhang Q, Fragnito N, Myers A, Sharma N. Plantarflexion moment prediction during the walking stance phase with an semg-ultrasound imaging-driven model. In: 2021 43rd Annual International Conference of the IEEE Engineering in Medicine & Biology Society (EMBC), 2021;6267–6272. IEEE

  45. Delp SL, Anderson FC, Arnold AS, Loan P, Habib A, John CT, Guendelman E, Thelen DG. Opensim: open-source software to create and analyze dynamic simulations of movement. IEEE Trans Biomed Eng. 2007;54(11):1940–50.

    PubMed  Article  Google Scholar 

  46. Kainz H, Hoang HX, Stockton C, Boyd RR, Lloyd DG, Carty CP. Accuracy and reliability of marker-based approaches to scale the pelvis, thigh, and shank segments in musculoskeletal models. J Appl Biomech. 2017;33(5):354–60.

    PubMed  Article  Google Scholar 

  47. Scott SH, Winter DA. A comparison of three muscle pennation assumptions and their effect on isometric and isotonic force. J Biomech. 1991;24(2):163–7.

    CAS  PubMed  Article  Google Scholar 

  48. Epstein M. Theoretical models of skeletal muscle: biological and mathematical considerations. New Jersey: Wiley; 1998. p. 52–3.

    Google Scholar 

  49. Guimaraes A, Herzog W, Hulliger M, Zhang Y, Day S. Effects of muscle length on the emg-force relationship of the cat soleus muscle studied using non-periodic stimulation of ventral root filaments. J Exp Biol. 1994;193(1):49–64.

    CAS  PubMed  Article  Google Scholar 

  50. Huijing PA. Important experimental factors for skeletal muscle modelling: non-linear changes of muscle length force characteristics as a function of degree of activity. Eur J Morphol. 1996;34(1):47–54.

    CAS  PubMed  Article  Google Scholar 

  51. Zajac FE. Muscle and tendon: properties, models, scaling, and application to biomechanics and motor control. Crit Rev Biomed Eng. 1989;17(4):359–411.

    CAS  PubMed  Google Scholar 

  52. Ao D, Song R, Gao J. Movement performance of human-robot cooperation control based on EMG-driven hill-type and proportional models for an ankle power-assist exoskeleton robot. IEEE Trans Neural Syst Rehabil Eng. 2017;25(8):1125–34.

    PubMed  Article  Google Scholar 

  53. Charlton IW, Tate P, Smyth P, Roren L. Repeatability of an optimised lower body model. Gait Posture. 2004;20(2):213–21.

    CAS  PubMed  Article  Google Scholar 

  54. Campanini I, Merlo A, Degola P, Merletti R, Vezzosi G, Farina D. Effect of electrode location on emg signal envelope in leg muscles during gait. J Electromyogr Kinesiol. 2007;17(4):515–26.

    CAS  PubMed  Article  Google Scholar 

  55. Farris DJ, Lichtwark GA. UltraTrack: Software for semi-automated tracking of muscle fascicles in sequences of B-mode ultrasound images. Comput Methods Programs Biomed. 2016;128:111–8.

    PubMed  Article  Google Scholar 

  56. Yamaguchi G, Sawa A, Moran D, Fessler M, Winters J. A survey of human musculotendon actuator parameters. In: Winters JM, Woo SL-Y, editors. Multiple muscle systems: biomechanics and movement organization. New York: Springer;1990. p. 717–73.

  57. MacLean MK, Ferris DP. Human muscle activity and lower limb biomechanics of overground walking at varying levels of simulated reduced gravity and gait speeds. PloS One. 2021;16(7):0253467.

    Article  CAS  Google Scholar 

  58. Liu MM, Herzog W, Savelberg HHCM. Dynamic muscle force predictions from EMG: an artificial neural network approach. J Electromyogr Kinesiol. 1999;9(6):391–400.

    CAS  PubMed  Article  Google Scholar 

  59. Barbero M, Merletti R, Rainoldi A. Atlas of muscle innervation zones: understanding surface electromyography and its applications. Berlin: Springer; 2012.

    Book  Google Scholar 

  60. Aagaard P, Andersen JL, Dyhre-Poulsen P, Leffers A-M, Wagner A, Magnusson SP, Halkjær-Kristensen J, Simonsen EB. A mechanism for increased contractile strength of human pennate muscle in response to strength training: changes in muscle architecture. J Physiol. 2001;534(2):613–23.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  61. Schmitz A, Silder A, Heiderscheit B, Mahoney J, Thelen DG. Differences in lower-extremity muscular activation during walking between healthy older and young adults. J Electromyogr Kinesiol. 2009;19(6):1085–91.

    PubMed  Article  Google Scholar 

  62. Mademli L, Arampatzis A, Morey-Klapsing G, Brüggemann G-P. Effect of ankle joint position and electrode placement on the estimation of the antagonistic moment during maximal plantarflexion. J Electromyogr Kinesiol. 2004;14(5):591–7.

    PubMed  Article  Google Scholar 

Download references

Acknowledgements

The authors would like to thank Alison Mayers and Mayank Singh for their help in experiments’ conduction and data collection. We would also like to thank all participants for their involvement in the study.

Funding

This research was funded by NSF Career Award Grant Number 2002261.

Author information

Authors and Affiliations

Authors

Contributions

QZ, JRF, and NS conceived and designed the experiments. QZ led the experimental setup, data collection, and data processing, analyzed and interpreted data, generated figures, and wrote the manuscript. NF collaborated with the experimental setup, processed ultrasound imaging data, and revised the manuscript. JRF provided technical support of using UltraTrack to process ultrasound imaging data and revised the manuscript. NS was the principal investigator, led the project, interpreted data, and revised the manuscript. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Nitin Sharma.

Ethics declarations

Ethics approval and consent to participate

The study was conducted according to the guidelines of the Declaration of Helsinki, and all experiments were approved by the Institutional Review Board of North Carolina State University (protocol number 20602 and date of approval 02/10/2020). Participants provided written and informed consent.

Consent for publication

Participants provided written and informed consent for publication.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Table S1.

Participants’ anthropometric characteristics and peak net PF moment during treadmill walking at each speed. Table S2. HNMs prediction performance under single-speed modes and inter-speed mode on Sub01. (RMSE unit: Nm). Table S3. HNMs prediction performance under single-speed modes and inter-speed mode on Sub02. (RMSE unit: Nm). Table S4. HNMs prediction performance under single-speed modes and inter-speed mode on Sub03. (RMSE unit: Nm). Table S5. HNMs prediction performance under single-speed modes and inter-speed mode on Sub04. (RMSE unit: Nm). Table S6. HNMs prediction performance under single-speed modes and inter-speed mode on Sub05. (RMSE unit: Nm). Figure S1. Right ankle joint kinematics and kinetics results from Participant Sub01. Figure S2. Right ankle joint kinematics and kinetics results from Participant Sub02. Figure S3. Right ankle joint kinematics and kinetics results from Participant Sub03. Figure S4. Right ankle joint kinematics and kinetics results from Participant Sub04. Figure S5. Right ankle joint kinematics and kinetics results from Participant Sub05. Figure S6. sEMG raw data and processed data from both LGS and SOL muscles during 20 seconds walking experiments on Participant Sub01. Figure S7. sEMG raw data and processed data from both LGS and SOL muscles during 20 seconds walking experiments on Participant Sub02. Figure S8. sEMG raw data and processed data from both LGS and SOL muscles during 20 seconds walking experiments on Participant Sub03. Figure S9. sEMG raw data and processed data from both LGS and SOL muscles during 20 seconds walking experiments on Participant Sub04. Figure S10. sEMG raw data and processed data from both LGS and SOL muscles during 20 seconds walking experiments on Participant Sub05. Figure S11. US imaging-derived MT and processed sEMG signals during the walking stance phase from Sub01. Figure S12. US imaging-derived MT and processed sEMG signals during the walking stance phase from Sub02. Figure S13. US imaging-derived MT and processed sEMG signals during the walking stance phase from Sub03. Figure S14. US imaging-derived MT and processed sEMG signals during the walking stance phase from Sub04. Figure S15. US imaging-derived MT and processed sEMG signals during the walking stance phase from Sub05. Figure S16. The Ultratrack toolbox interface with US imaging-derived MT tracking procedures and results without key-frame correction on Participant Sub07 walking at 1.00 m/s. Figure S17.The Ultratrack toolbox interface with US imaging-derived MT tracking procedures and results with key-frame correction on Participant Sub07 walking at 1.00 m/s. Figure S18. Exampled calibration and prediction performance on Participant Sub01 by using five HNM categories. The blue, red, orange, purple, and green bars represent the results by using the sEMG-, US imaging-, sEMG-US imaging-, random signal (RS)-, and sEMG-RS-driven HNMs under both single-speed modes and inter-speed mode. (a) Calibration RMSE and prediction RMSE values under different speed modes. (b) Calibration R-square and prediction R-square values under different speed modes. Figure S19. Exampled results of the ankle joint net moment and processed sEMG signals from LGS, MGS, SOL, and TA muscles throughout the recorded walking duration. Data are collected during the treadmill walking experiments on Sub01, where muscle activation levels from LGS and MGS muscles are comparable, and the activation level of the TA muscle is relatively small during the walking stance phase. (a) Walking speed at 0.5 m/s. (b) Walking speed at 0.75 m/s.

Additional file 2. The LGS and SOL muscles' thickness tracking results from sequential US images on Participant Sub01 during treadmill walking at 0.50 m/s.

Additional file 3. The LGS and SOL muscles' thickness tracking results from sequential US images on Participant Sub01 during treadmill walking at 0.75 m/s.

Additional file 4. The LGS and SOL muscles' thickness tracking results from sequential US images on Participant Sub01 during treadmill walking at 1.00 m/s.

Additional file 5. The LGS and SOL muscles' thickness tracking results from sequential US images on Participant Sub01 during treadmill walking at 1.25 m/s.

Additional file 6. The LGS and SOL muscles' thickness tracking results from sequential US images on Participant Sub01 during treadmill walking at 1.50 m/s.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Zhang, Q., Fragnito, N., Franz, J.R. et al. Fused ultrasound and electromyography-driven neuromuscular model to improve plantarflexion moment prediction across walking speeds. J NeuroEngineering Rehabil 19, 86 (2022). https://doi.org/10.1186/s12984-022-01061-z

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12984-022-01061-z

Keywords

  • Human intent prediction
  • Neuromuscular model
  • B-mode ultrasound imaging
  • Surface electromyography
  • Rehabilitative/assistive robotics
  • Sensor fusion