Skip to main content

ORIGINAL RESEARCH article

Front. Bioeng. Biotechnol., 07 September 2022
Sec. Biomechanics

EMG-driven musculoskeletal model calibration with estimation of unmeasured muscle excitations via synergy extrapolation

  • 1Rice Computational Neuromechanics Lab, Department of Mechanical Engineering, Rice University, Houston, TX, United States
  • 2Biomechanics, Rehabilitation, and Integrative Neuroscience (BRaIN) Lab, VA Northern California Health Care System, Martinez, CA, United States
  • 3Department of Physical Medicine and Rehabilitation, Davis School of Medicine, University of California, Sacramento, CA, United States

Subject-specific electromyography (EMG)-driven musculoskeletal models that predict muscle forces have the potential to enhance our knowledge of internal biomechanics and neural control of normal and pathological movements. However, technical gaps in experimental EMG measurement, such as inaccessibility of deep muscles using surface electrodes or an insufficient number of EMG channels, can cause difficulties in collecting EMG data from muscles that contribute substantially to joint moments, thereby hindering the ability of EMG-driven models to predict muscle forces and joint moments reliably. This study presents a novel computational approach to address the problem of a small number of missing EMG signals during EMG-driven model calibration. The approach (henceforth called “synergy extrapolation” or SynX) linearly combines time-varying synergy excitations extracted from measured muscle excitations to estimate 1) unmeasured muscle excitations and 2) residual muscle excitations added to measured muscle excitations. Time-invariant synergy vector weights defining the contribution of each measured synergy excitation to all unmeasured and residual muscle excitations were calibrated simultaneously with EMG-driven model parameters through a multi-objective optimization. The cost function was formulated as a trade-off between minimizing joint moment tracking errors and minimizing unmeasured and residual muscle activation magnitudes. We developed and evaluated the approach by treating a measured fine wire EMG signal (iliopsoas) as though it were “unmeasured” for walking datasets collected from two individuals post-stroke–one high functioning and one low functioning. How well unmeasured muscle excitations and activations could be predicted with SynX was assessed quantitatively for different combinations of SynX methodological choices, including the number of synergies and categories of variability in unmeasured and residual synergy vector weights across trials. The two best methodological combinations were identified, one for analyzing experimental walking trials used for calibration and another for analyzing experimental walking trials not used for calibration or for predicting new walking motions computationally. Both methodological combinations consistently provided reliable and efficient estimates of unmeasured muscle excitations and activations, muscle forces, and joint moments across both subjects. This approach broadens the possibilities for EMG-driven calibration of muscle-tendon properties in personalized neuromusculoskeletal models and may eventually contribute to the design of personalized treatments for mobility impairments.

1 Introduction

Muscle force is an important biomechanical variable for both research and clinical purposes. Knowledge of muscle forces during a variety of motion tasks could facilitate the development of more effective treatments for neuromusculoskeletal disorders such as stroke (Shao et al., 2009), Parkinson’s disease (Cano-de-la-Cuerda et al., 2010), and knee osteoarthritis (Kim et al., 2009). Specifically, such knowledge could lead to an improved understanding of the control strategies employed by the central neural systems (CNS) (Contessa and Luca, 2013; Del Vecchio et al., 2018), internal biomechanical quantities such as joint contact forces (Correa et al., 2010; Sasaki and Neptune, 2010; Manal and Buchanan, 2013; Walter et al., 2014; Serrancolí et al., 2016; Hoang et al., 2019), and external biomechanical quantities such as joint moments generated by muscles (Manal and Buchanan, 2003; Buchanan et al., 2005; Sartori et al., 2012; Meyer et al., 2017). While researchers are seeking to develop new experimental methods to measure muscle or tendon forces in vivo during human movement (e.g., Martin et al., 2018), direct measurement of muscle forces in vivo remains inherently challenging, costly, and ethically problematic to perform.

This situation has motivated the development of computational approaches that can estimate muscle forces from subject movement data. The two primary computational approaches used for this purpose are static optimization (SO) and EMG-driven musculoskeletal modeling. Both approaches typically utilize a geometric musculoskeletal model actuated by Hill-type muscle-tendon models (Zajac, 1989), where the control inputs to the muscle-tendon models can be either muscle excitations (equivalent to processed experimental-electromyographic (EMG) data) or muscle activations (muscle excitations after being passed through an activation dynamics model). SO is used to estimate muscle activations when EMG data are wholly (Crowninshield and Brand, 1981; Anderson and Pandy, 2001; Heintz and Gutierrez-Farewik, 2007) or partially (Sartori et al., 2014; Zonnino and Sergi, 2019) missing from important modeled muscles. This approach resolves the muscle redundancy problem by using nonlinear optimization to adjust the predicted muscle activations such that the sum of squares (or some other power) of muscle activations is minimized and the predicted net joint moments match the inverse dynamic net joint moments (Anderson and Pandy, 2001; Ackermann and van den Bogert, 2010). In contrast, EMG-driven modeling is used to estimate muscle excitations when EMG data are wholly or mostly available from important modeled muscles (Lloyd and Besier, 2003; Manal and Buchanan, 2003; Amarantini and Martin, 2004; Shao et al., 2009; Menegaldo et al., 2014; Sartori et al., 2014; Pizzolato et al., 2015; Meyer et al., 2017). To resolve the muscle redundancy problem, this approach uses processed EMG data to define the shapes of the predicted muscle excitations and then uses nonlinear optimization to adjust activation dynamics and Hill-type muscle-tendon model parameter values such that the sum of squares of errors between predicted and inverse dynamic net joint moments is minimized. While SO is extremely fast computationally, it can underestimate muscle activations since co-activation between agonist and antagonist muscles is minimized (Herzog and Binding, 1992), and it can produce unrealistic abrupt activation changes since the optimization process solves each time frame independently (Vilimek, 2007; Schellenberg et al., 2015). More importantly, it does not provide a way to calibrate activation dynamics and muscle-tendon model parameter values to a subject’s movement data, which can adversely affect the reliability of the estimated muscle forces (Serrancolí et al., 2016).

While calibration of activation dynamics and muscle-tendon model parameter values is built into the EMG-driven modeling process, several practical challenges exist with collecting EMG data from all muscles that contribute significantly to a specified movement task (Sartori et al., 2014; Péter et al., 2019; Zonnino and Sergi, 2019; Ao et al., 2020; Gurchiek et al., 2020). First, surface electrodes, which are non-invasive and easily applied, cannot measure EMG signals from deep muscles that contribute significantly to joint moments. Common examples are the iliacus and psoas muscles (Sartori et al., 2014), which significantly contribute to the hip flexion moment during walking. Second, though fine wire electrodes can measure EMG signals from deep muscles, they are invasive, require special skill and significant preparation time to insert, and may cause the subject discomfort and pain, limiting their utilization. Third, in some situations, deep muscles may not be reachable by any type of electrode. For example, use of a fine wire electrode may be contraindicated for safety reasons in subjects who have a cancerous tumor near an important deep muscle. Fourth, the number of available EMG channels is often less than that needed to drive a neuromusculoskeletal model for multi-joint movements (e.g., typically more than 10 channels per leg for walking and running). However, human movement labs often have only eight or 16 channels available for EMG recording. These challenges are important since missing EMG data from critical muscles may have a domino effect on the reliability of force estimates for other muscles that span the same joints (Pizzolato et al., 2015; Zonnino and Sergi, 2019).

Given the challenges described above, researchers have sought to develop various computational methods for estimating missing EMG signals during the EMG-driven model calibration process. One such method utilizes Gaussian process regression models to describe the synergistic relationship between a subset of muscles, enabling the estimation of unmeasured muscle excitations using information provided by a subset of measured muscle excitations (Gurchiek et al., 2020). However, muscle excitations associated with “unmeasured” muscles must be known initially to perform the necessary model training process, making this method infeasible when the “unmeasured” muscle excitations are truly unmeasurable due to experimental limitations or safety considerations. Another method utilizes low-dimension a sets of impulsive excitation primitives to estimate unmeasured muscle excitations (Neptune et al., 2009; Sartori et al., 2013; Pizzolato et al., 2015). Once excitation primitives are derived from measured muscle excitations, each muscle is assigned to a module by assessing associated weighting factors, where muscles without EMG measurements are assumed to belong to the same module as measured muscles that share the same innervation and contribute to the same mechanical action. Pizzolato et al. (2015) also minimally adjusted primitive-driven excitaitons for muscles with experimental EMG data to improve joint moment estimation in EMG-assisted mode. However, these adjustments masked the omission of active force generating properties for muscles without EMG data (i.e., iliacus and psoas), resulting in noticible hip joint moment prediction errors. Furthermore, none of these studies evaluated the accuracy of predicted umeasurd muscle excitations due to the lack of corresponding experimental EMG data.

More recently, the muscle synergy concept has been investigated for estimating muscle activations via SO or muscle excitations via EMG-driven modeling (Bianco et al., 2018; Ao et al., 2020; Michaud et al., 2020; Shourijeh and Fregly, 2020). A muscle synergy consists of a time-varying synergy excitation (or activation) and a corresponding time-invariant synergy vector containing weights that define how the synergy excitation (or activation) contributes to the excitation (or activation) of all muscles (Cappellini et al., 2006; Tresch et al., 2006; Ting and Chvatal, 2010; Banks et al., 2017). Muscle synergies are useful because they allow a large number of measured or modeled muscle excitations (or activations) to be represented by a small number of muscle synergies (typically between 3 and 6) (Tresch et al., 2006; Ivanenko et al., 2005; Ting and Chvatal, 2010; Banks et al., 2017). Michaud et al. (2020) and Shourijeh and Fregly (2020) imposed a synergy structure on muscle activations estimated via SO, which required solving for muscle activations over all time frames simultaneously. In both studies, muscle-tendon model parameter values were not calibrated simultaneously, no experimental EMG data were used to inform the muscle activation solutions, and the accuracy of predicted muscle activations compared to experimental EMG measurements were no better than from SO. In another recent study, Ao et al. (2020) used synergy excitations calculated from measured muscle excitations to predict synergy vector weights for unmeasured muscle excitations via a simplied EMG-driven modeling process. To evaluate the feasibility of the approach (which was termed “synergy extrapolation” or SynX), the authors used an EMG-driven model whose activation dynamics and muscle-tendon model parameters were already calibrated to subject movement data using a complete set of EMG measurements. Predictions of muscle excitations that were measured experimentally using fine wire electrodes but treated as unmeasured for SynX evaluation purposes showed excellent agreement with corresponding EMG measurements. However, since a pre-calibrated EMG-driven model was used for the evaluation, it remains unknown whether SynX can predict unmeasured muscle excitations reliably when EMG-driven model calibration is performed simultaneously.

This study extends the capabilities of SynX so that it can estimate missing EMG signals while simultaneously calibrating activation dynamics and muscle-tendon model parameter values in an EMG-driven model. The approach was developed and evaluated using gait datasets collected from two subjects post-stroke performing treadmill walking at self-selected and fastest-comfortable speeds. EMG signals measured bilaterally from iliopsoas using fine-wire electrodes were treated as “unmeasured” and used to evaluate the reliability of the method quantitatively. The computational approach uses nonlinear optimization to calibrate three categories of design variables simultaneously: 1) EMG-driven model parameters, 2) synergy vector weights and average values for constructing unmeasured excitations for muscles without associated EMG data, and 3) synergy vector weights and average values for constructing residual excitations for muscles with associated EMG data. The cost function was formulated as a trade-off between joint moment matching accuracy and unmeasured and residual muscle activation minimization. Different methodological choices, including the number of synergies and variability of synergy vector weights across trials, were investigated to determine the best choices for analyzing experimentally measured walking motions and generating computationally predicted walking motions. EMG-driven lower extremity models were calibrated for both subjects using the standard method with no missing EMG signals and SynX with missing iliacus and psoas EMG signals. Muscle excitations, activations, forces, and net joint moments along with EMG-driven model parameter values produced by the two calibration methods were compared for the walking trials used in the calibration process and additional walking trials held back for validation purposes.

2 Methods

2.1 Experimental data collection and processing

Previously published walking datasets collected from a high-functioning hemiparetic subject (S1, male, height 1.70 m, mass 80.5 kg, age 79 years, right-hand hemiparesis, lower extremity Fugl-Meyer Motor Assessment score of 32 out of a maximum 34) and a low-functioning hemiparetic subject post-stroke (S2, male, height 1.83 m, mass 88.5 kg, age 62 years, right-hand hemiparesis, lower extremity Fugl-Meyer Motor Assessment score of 25 out of a maximum 34) were used for this study (Meyer et al., 2017; Li et al., 2020). Video motion capture (Vicon Corp., Oxford, United Kingdom), ground reaction (Bertec Corp., Columbus, OH, United States), and EMG (Motion Lab Systems, Baton Rouge, LA, United States) data were recorded simultaneously while subjects walked on a split-belt instrumented treadmill (Bertec Corp., Columbus, OH, United States) at their self-selected (0.5 m/s for S1 and 0.45 m/s for S2) and fastest-comfortable (0.8 m/s for S1 and 0.65 m/s for S2) speeds. EMG data consisted of sixteen channels collected from each leg at 1,000 Hz using a combination of surface and fine wire electrodes (Table 1), which made EMG data available for important deep muscles such as iliopsoas. All experimental procedures were approved by the University of Florida Health Science Center Institutional Review Board (IRB-01), and both subjects provided written informed consent before participation. Raw motion capture and ground reaction data were low-pass filtered with a cut-off frequency of 7/tf Hz, where tf is the period of the gait cycle, while raw EMG data were high-pass filtered at 40 Hz, demeaned, full-wave rectified, low-pass filtered at 3.5/tf Hz, and normalized to maximum values over all experimental gait cycles (McLean et al., 2005; Meyer et al., 2017). Henceforth processed EMG data will be referred to as “muscle excitations.” Thirty-four and thirty-three muscles in each leg of the model that had either surface or fine-wire EMG data were kept for analysis for S1 and S2, respectively (Table 1). Data from ten gait cycles (five cycles per speed) per leg were randomly selected for EMG-driven model calibration. After pre-processing, data from each gait cycle were resampled to 101-time points representing initial heel-strike (0%) to subsequent heel-strike (100%) of the same foot. Twenty additional time frames before the start of each gait cycle were retained to account for a maximum electromechanical delay of approximately 100 ms for each muscle, which made each gait cycle possess 121 time points.

TABLE 1
www.frontiersin.org

TABLE 1. List of muscles in the model for each subject and which DOF each muscle actuates.

2.2 Musculoskeletal model analyses

In preparation for EMG-driven model calibration, we performed a sequence of five musculoskeletal model analyses. First, a generic full-body OpenSim musculoskeletal model (Rajagopal et al., 2016) was scaled to match each subject’s anthropometry using OpenSim 4.0 (Delp et al., 2007; Seth et al., 2018). Each leg of the model possessed six degrees of freedom (DOFs), including hip flexion/extension (HipFE), hip adduction/abduction (HipAA), hip internal/external rotation (HipRot), knee flexion/extension (KneeFE), ankle plantarflexion/dorsiflexion (AnklePD), and ankle inversion/eversion (AnkleIE). Second, the locations and orientations of joint centers and axes, respectively, for the hip, knee, and ankle of each leg were adjusted via nonlinear optimization such that surface markers on the OpenSim model tracked experimentally measured surface marker positions as closely as possible for isolated joint motion and walking trials (Reinbolt et al., 2005). Third, inverse kinematic (IK) analyses were performed with OpenSim using experimental marker motion data from the walking trials to obtain the time histories of lower body joint angles. Fourth, inverse dynamic (ID) analyses were performed with OpenSim using IK joint angles and experimental ground reaction data for the same walking trials to calculate lower extremity joint moments. Fifth, for each muscle-tendon actuator in the model, a set of surrogate musculoskeletal geometric models was fitted to approximate the subject’s muscle-tendon lengths, velocities, and moment arms as a function of lower extremity joint angles and angular velocities (Menegaldo et al., 2004; Meyer et al., 2017).

2.3 EMG-driven model calibration with SynX

The proposed computational method employing synergy extrapolation (SynX) to estimate missing muscle excitations during EMG-driven model calibration involves five steps (Figure 1), where all steps except the first one are performed within a nonlinear optimization process. First, muscle synergy analysis (MSA) is performed on measured muscle excitations via principal component analysis (PCA) to extract time-varying synergy excitations (henceforth referred to as “measured synergy excitations”) along with time-invariant synergy vector weights that define how each measured synergy excitation contributes to all measured muscle excitations (see details in Extraction of measured synergy excitations). Second, unmeasured muscle excitations are constructed by linearly combining the measured synergy excitations using the current guesses for the unmeasured synergy vector weights and average values (see details in Estimation of unmeasured muscle excitations). Third, residual muscle excitations added to the experimental muscle excitations are constructed by linearly combining measured synergy excitations using the current guesses for the residual synergy vector weights and average values (see details in Estimation of residual muscle excitations). Fourth, muscle activations, forces, and net joint moments are calculated by the EMG-driven musculoskeletal modeling using experimental joint kinematics and the current guesses for muscle excitations and muscle-tendon model parameter values when residual excitations are and are not included in the calculation of net joint moments (see details in Formulation of EMG-driven musculoskeletal model). Fifth, the nonlinear optimization adjusts all synergy and muscle-tendon model parameters to reduce the multi-objective cost function (see details in Calibration of EMG-driven musculoskeletal model). Below we describe each of these five steps as applied to 10 cycles of gait data for each leg of the two experimental subjects.

FIGURE 1
www.frontiersin.org

FIGURE 1. Flow chart of EMG-driven model calibration with synergy extrapolation (SynX) to estimating missing EMG signals. Orange rectangle boxes indicate optimization design variables, including activation dynamics model parameters, muscle-tendon model parameters, average values and synergy vector weights for unmeasured muscle excitations, and average values and synergy vector weights for residual excitations. The design variables were calibrated simultaneously by solving a nonlinear multi-objective optimization problem, described as “Params + SynX + Res” calibration. The objective function was formulated as a trade-off between joint moment tracking accuracy (with and without residual excitations included) and the magnitude of estimated unmeasured and residual muscle activations. Green parallelogram boxes indicate variables that were measured experimentally. Blue rectangle boxes indicate important intermediate variables generated within the SynX process. MSA stands for muscle synergy analysis.

2.3.1 Extraction of measured synergy excitations

Measured synergy excitations were extracted from the measured muscle excitations (excluding iliopsoas excitations when using SynX) by performing muscle synergy analysis (MSA). A previous study explored how methodological choices involved in MSA affect synergy extrapolation performance when a pre-calibrated EMG-driven musculoskeletal model is used (Ao et al., 2020). That study demonstrated that principal component analysis (PCA) with five or six synergies consistently predicted unmeasured muscle excitations accurately, while non-negative matrix factorization (NMF) could not achieve acceptable prediction accuracy. Moreover, EMG normalization did not significantly affect synergy extrapolation performance. Thus, in this study, measured muscle excitations were normalized to their maximum values over all trials, and PCA was selected for performing MSA. For between four and seven synergies, measured muscle excitations were decomposed via PCA and represented as:

em= WmHm+μm+εm(1)

where em stands for measured muscle excitations, Wm for time-varying measured synergy excitations, and Hm for associated measured synergy vector weights. εm accounts for the decomposition residuals that could not be explained by WmHm+μm, where μm specifies the average value of each measured muscle excitation.

MSA was performed to extract measured synergy excitations based on three assumptions about how measured synergy vector weights vary across trials. First, measured synergy vector weights were assumed to be trial-specific, necessitating a separate PCA decomposition for each trial individually, where em was a ntrial (frames per trial) × k (measured excitations) matrix. Wm became an ntrial (frames per trial) × p (synergies) matrix and Hm became a = p (synergies) × k (measured excitations) matrix. Second, measured synergy vector weights were assumed to be speed-specific, where the measured muscle excitations from trials with the same speed were concatenated such that em became nspeed (frames per speed) × k (measured excitations). Thus, Wm became nspeed (frames per speed) × p (synergies) and Hm became p (synergies) × k (measured excitations). Third, measured synergy vector weights were assumed to be subject-specific. For this assumption, PCA was run on a concatenated matrix of measured muscle excitations over all trials, where  em was nsubject (frames over all trials) by k (measured excitations). Wm had dimension of nsubject (frames over all trials) × p (synergies) and Hm had the dimension of p (synergies) × k (measured excitations). In this study, five trials were analyzed for each of the two walking speeds, which led nspeed to be 605 [121 (frames per trial) × 5 (trials per speed)] and nsubject to be 1,210 [121 (frames per trial) × 10 (trials over two speeds)]. PCA-based MSAs were performed using the “pca” command in MATLAB (The Mathworks, Natick, MA).

2.3.2 Estimation of unmeasured muscle excitations

Following MSA, unmeasured muscle excitations were constructed from the measured synergy excitations Wm as shown in Eq. 2 below:

eSynX= WmHSynX+μSynX(2)

where eSynX represents the unmeasured muscle excitations, HSynX represents the unmeasured synergy vector weights, and μSynX represents the average value of each unmeasured muscle excitation. Both HSynX and μSynX were design variables in the optimization problems formulated within our EMG-driven calibration approach (see more details in Section 2.3.5). Unmeasured synergy vector weights were assumed to be trial-specific, speed-specific, or subject-specific (henceforth referred to as “categories of unmeasured synergy vector weights”), where Wm  used for reconstruction was consistent with the corresponding category. More specifically, ten sets of HSynX and μSynX were generated for the trial-specific category, two sets of HSynX  and μSynX for the speed-specific category, and only one set of HSynX and μSynX for the subject-specific category. For all categories of unmeasured synergy vector weights, HSynX was of dimensions p (synergies) × q (unmeasured excitations) and μSynX was 1 × q (unmeasured excitations).

2.3.3 Estimation of residual muscle excitations

Similar to unmeasured muscle excitations, residual muscle excitations were constructed from the measured synergy excitations Wm using the following relationship:

eres= WmHres+μres(3)

where eres denotes residual muscle excitations, Hres denotes residual synergy vector weights, and μres denotes the average value of the residual muscle excitation. Again, residual synergy vector weights were assumed to be trial-specific, speed-specific, or subject-specific specific (henceforth referred to as “categories of residual synergy vector weights”), where Wm  used for reconstruction was consistent with the corresponding category. For all categories of residual synergy vector weights, Hres was of dimensions p (synergies) × q (measured excitations) and μres was 1 × q (measured excitations). The estimated residual muscle excitations were added to the measured muscle excitations to produce adjusted measured excitations:

emres= em+eres(4)

where emres refers to measured muscle excitations with residual excitations included.

2.3.4 Formulation of EMG-driven musculoskeletal model

Once unmeasured and residual muscle excitations were constructed, an EMG-driven musculoskeletal model was used to predict the six lower extremity net joint moments (Lloyd and Besier, 2003; Manal and Buchanan, 2003; Amarantini and Martin, 2004; Shao et al., 2009; Menegaldo et al., 2014; Sartori et al., 2014; Pizzolato et al., 2015; Meyer et al., 2017; Ao et al., 2020). First, muscle excitations e(t) were scaled by muscle-specific scale factors sEMG between 0.05 and 1 to reflect the fact that true maximum excitation levels are likely to be higher than those observed experimentally. Next, neural activation u(t) was calculated from muscle excitation e(t) using a published activation dynamics model (He et al., 1991) that used a first-order ordinary differential equation to define the relationship between e(t) and u(t):

du(t)dt=(c1e(td)+c2)(e(td)u(t))(5)
c1=1/τact1/τdact (6)
c2=1/τdact (7)

where τact and τdact are activation and deactivation time constants, respectively, τdact was assumed to be 4τact (Zajac, 1989; Meyer et al., 2017), and d denotes an electromechanical time delay. Once neural activation was computed, a nonlinear function was used to determine the corresponding muscle activation a(t) (Manal and Buchanan, 2003):

a(t)=(1c3)u(t)+c3[g1g2(u(t)+g3)g4+g5+1](8)

where c3 is an activation nonlinearity constant, g1 to g5 are constant coefficients determined by fitting published experimental data from isometric contractions (Manal and Buchanan, 2003).

Two sets of muscle activations were calculated when residual muscle excitations were and were not included, respectively. The muscle activations predicted from a combination of unmeasured muscle excitations (eSynX) and measured muscle excitations (em) were defined as apre, while the muscle activations predicted from a combination of unmeasured muscle excitations (eSynX) and measured muscle excitations with residual excitations included (emres) were defined as apreres. Residual muscle activations ares were then defined as the difference between apreres and apre:

ares= apreres apre(9)

Next, taking both apre and apreres as inputs, a Hill-type muscle-tendon model with rigid tendon (Hill, 1938; Zajac, 1989; Meyer et al., 2017) was used to predict the force generated by a given muscle-tendon actuator as shown below:

F(t)=FoM·[a(t)·fl(l˜M(t))·fv(v˜M(t))+fp(l˜M(t))]cosα(10)

where F(t) is the force generated by the muscle-tendon actuator, FoM is the maximum isometric force of the muscle, a is the muscle activation, l˜M(t) and v˜M(t) are the time-varying normalized muscle fiber length and velocity, respectively, and α is the pennation angle of the muscle. fl(l˜M(t)) and fv(v˜M(t)) define the normalized active muscle force-length and force-velocity relationships, respectively, while fp(l˜M(t)) defines the normalized passive muscle force-length relationship. l˜M(t) and v˜M(t) were calculated using the following equations, which assume a rigid tendon:

l˜M(t)=lMT(t)lsTloM(11)
v˜M(t)=vMT(t)10·loM(12)

where lMT is muscle-tendon length and vMT muscle-tendon velocity, loM is optimal muscle fiber length, and lsT is tendon slack length. With this muscle-tendon model, muscle forces were estimated when residual excitations were included (Fpreres) and not included (Fpre) during prediction respectively.

Once Fpreres and Fpre were computed for all muscles in the model, their contributions to net joint moments M were calculated using the corresponding muscle moment arms:

M(t)=F(t)· r(t)(13)

where r is muscle moment arm defined as the negative of the partial derivative of muscle-tendon length lMT with respect to generalized coordinate θ (An et al., 1984):

r(t)=lMT(t)θ(14)

The negative sign in Eq. 14 was implemented for consistency with the OpenSim modeling environment. As required by the cost function for EMG-driven model calibration, net joint moments were calculated when residual excitations were (Mpreres) and were not (Mpre) included in the measured muscle excitations.

2.3.5 Calibration of EMG-driven musculoskeletal model

Calibration of the EMG-driven musculoskeletal model with simultaneous estimation of missing iliacus and psoas EMG signals was performed by using nonlinear optimization to adjust four categories of design variables: 1) activation dynamics model parameters consisting of EMG scale factor, sEMG, electromechanical delay d, activation time constant τact, and activation nonlinearity constant c3, 2) scaling factors for muscle-tendon model parameters consisting of optimal muscle fiber length loM and tendon slack length lsT, 3) synergy vector weights HSynX and average values μSynX associated with unmeasured muscle excitations, and 4) synergy vector weights Hres and average values μres associated with residual muscle excitations. To develop and evaluate the performance of our SynX approach, we formulated three optimization problems using various combinations of design variables (Table 2). First, to estimate missing EMG signals, we calibrated EMG-driven model parameters, synergy vector weights plus average values for unmeasured muscle excitations, and synergy vector weights plus average values for residual muscle excitations (termed “Params + SynX + Res,” see details in Case 1: Params + SynX + Res). Second, to explore how SynX performance is influenced by the inclusion of residual muscle excitations, we calibrated only EMG-driven model parameters and synergy vectors weights plus average values for unmeasured muscle excitations (termed “Params + SynX,” see details in Case 2: Params + SynX). Third, to evaluate how well important biomechanical variables (e.g., muscle activations, muscle forces, and net joint moments) can be estimated using “Params + SynX + Res” calibration, we calibrated only EMG-driven model parameters using the complete set of EMG signals with no muscle excitations predicted by SynX (termed “Params,” see details in Case 3: Params). For all three optimization problems, if unmeasured or residual muscle excitations were needed, we explored all possible methodological combinations of the number of synergies, category of unmeasured synergy vector weights, and category of residual synergy vector weights.

TABLE 2
www.frontiersin.org

TABLE 2. Summary of calibration cases, which were named based on categories of design variables included in the optimization problem formulation.

2.3.5.1 Case 1: Params + SynX + Res

EMG-driven model calibration typically adjusts muscle forces by altering muscle-tendon model parameter values such that the difference between model-predicted and inverse dynamic (ID) joint moments are minimized. However, when unmeasured muscle excitations are estimated via SynX during EMG-driven model calibration, four terms are minimized simultaneously: 1) sum of squares of errors between model-predicted (Mpreres) and inverse dynamic (Mexp) joint moments when residual muscle excitations are included in joint moment calculations, termed JMomres; 2) sum of squares of errors between model-predicted (Mpre) and inverse dynamic (Mexp) joint moments when residual muscle excitations are not included in the joint moment calculations, termed JMom; 3) sum of squares of unmeasured muscle activations estimated by SynX (aSynX), termed JactSynX; and 4) sum of squares of residual muscle activations (ares), termed Jactres. Thus, the cost function for EMG-driven model calibration with SynX was formulated as:

minJJMomres+JMom +JactSynX +Jactres (15)

where

JMomres=(MpreresMexpMAD1 )2(16)
JMom=(MpreMexp MAD2)2(17)
JactSynX=(aSynXMAD3 )2(18)
Jactres=(aresMAD4 )2(19)

All four cost function terms (1619) were normalized by a maximum allowable deviation (MAD). A sensitivity analysis was performed to determine a MAD value for each cost function term, as described in the Appendix. Specific details about initial guesses and upper/lower bounds for each design variable, additional constraints, and cost function penalty terms can be found in Supplementary Table S1 in the Appendix and previously published papers (Meyer et al., 2017; Ao et al., 2020). All optimizations were performed using MATLAB’s “fmincon” function with the sequential quadratic programming algorithm.

2.3.5.2 Case 2: Params + SynX

To assess how well unmeasured muscle excitations can be predicted when no residual muscle excitations are included, we estimated unmeasured muscle excitations as in Case 1 but without adding residual muscle excitations to the measured muscle excitations. For this case, the cost function was:

minJ(MpreMexp MAD2)2 +(aSynXMAD3 )2(20)

where Mpre denotes the joint moments calculated from a combination of measured (em) and unmeasured (eSynX) muscle excitations.

2.3.5.3 Case 3: Params

To evaluate the performance of our computational approach for different methodological choices, we performed EMG-driven model calibration using the full set of EMG signals, where no EMG signals were treated as “unmeasured,” and only activation dynamics and muscle-tendon model parameter values were optimized to match the experimental joint moments from inverse dynamics:

minJ(MpreMexpMAD2 )2(21)

where Mpre represents model-predicted net joint moments produced by a complete set of muscle excitations.

2.4 EMG-driven model evaluation

Several common metrics were used to evaluate the outcome of muscle synergy analysis and the reliability of unmeasured muscle excitations generated using SynX with all possible methodological combinations. First, the variance accounted for (VAF) was computed between experimental and reconstructed muscle excitations of measured muscles, where the number of synergies in each leg was determined using a threshold criterion of 95% VAF (Tresch et al., 2006; Steele et al., 2013). Next, the EMG-driven model calibration process utilizing SynX was evaluated in three stages. In the first stage, the influence of different methodological choices (i.e., number of synergies, category of unmeasured and residual synergy vector weights) on SynX performance was investigated, and the choices that produced the most reliable estimates of unmeasured muscle excitations were identified. The two “most reliable” SynX methodological combinations for two distinct situations were identified by quantifying how well unmeasured muscle excitations and activations for psoas and iliacus could be predicted. Root mean square error (RMSE) and Pearson correlation coefficient r between experimental (from “Params” case) and estimated (from “Params + SynX + Res” case) muscle excitations and activations for iliacus and psoas across two speeds were computed to quantify matching of magnitude and shape, respectively. Correlation was interpreted quantitatively as weak (r < 0.35), moderate (0.35 < r ≤ 0.67), strong (0.67 < r ≤ 0.9), or very strong (r ≤ 0.9) (Taylor, 1990). The first situation involved predictions made when synergy vector weights and average values for unmeasured muscle excitations, along with muscle-tendon model parameter values, were calibrated via SynX. This situation (called “calibration”) is how the EMG-driven model would analyze walking trials used in the calibration process, where synergy vector weights and average values for unmeasured and residual muscle excitations can be calibrated on a trial-specific, speed-specific, or subject-specific basis. The second situation involved predictions made when synergy vector weights and average values for unmeasured muscle excitations, along with muscle-tendon model parameter values, were pre-calibrated via SynX. This situation (called “validation”) is how the EMG-driven model would be used to analyze experimental walking trials not used in the calibration process or to generate computationally predicted walking motions, where pre-calibrated synergy vector weights and average values for unmeasured muscle excitations must be utilized on either a speed-specific or subject-specific basis (Meyer et al., 2016; Sauder et al., 2019). In the second stage, muscle activations, forces, and net joint moments from “Params + SynX + Res” calibration were compared with those from “Params” calibration for both subjects using the first “most reliable” methodological combination applied to 10 calibration walking trials (five per speed). In the third stage, muscle activations, forces, and net joint moments from “Params + SynX + Res” validation were compared with those from “Params” validation for both subjects using the second “most reliable” methodological combination applied to 10 validation walking trials (five per speed) not utilized in the calibration process. For this stage, unmeasured synergy vectors weights and average values were fixed at their calibrated values. Mean absolute errors (MAE) between inverse dynamic and model-predicted net joint moments were calculated for “Params + SynX + Res” and “Params” cases using calibration trials in the second stage and validation trials in the third stage, respectively. Variance accounted for (VAF) was computed between experimental and model-predicted unmeasured muscle excitations, where the number of synergies in each leg was determined using a threshold criterion of 95% VAF (Tresch et al., 2006; Steele et al., 2013).

Multiple statistical analyses were performed to assess whether the calculated metrics resulting from different SynX methodological choices were statistically different. To assess whether methodological choices for MSA had a statistically significant impact on the reconstruction performance of measured muscle excitations, we performed a two-factor ANOVA with a Tukey-Kramer post-hoc analysis on VAF values. Second, to compare SynX performance among different SynX methodological choices, we performed a three-factor ANOVA tests on r and RMSE values across both patients and all calibration trials, followed by paired t-tests for comparing categories of synergy vector weights for a given number of synergies. For example, when we compared the difference among categories of unmeasured synergy vector weights, we paired r or RMSE values such that each pair shared the same category of residual synergy vector weights and number of synergies when associated with the same leg. Third, to investigate whether the inclusion of residual muscle excitations influenced SynX performance, we performed paired t-tests on r and RMSE values between “Params + SynX + Res” calibration with each category of residual synergy vector weights and “Params + SynX” calibration with no residual excitations. Fourth, we conducted paired t-tests to compare the joint moment matching errors (MAE values) between “Params” calibration and “Params + SynX + Res” calibration, where residual excitations were used in the calibration process but not included in the final joint moment calculation. All statistical analyses were performed in MATLAB, and the level of statistical significance was set at p < 0.05.

3 Results

3.1 Analysis of experimental muscle synergies

The two-way ANOVA applied to mean VAF values between reconstructed and experimental muscle excitations of measured muscles revealed the main effects of the number of synergies (p < 0.05) and the category of measured synergy vector weights (p < 0.05) on the variance explained by the factorization of measured muscle excitations. Post-hoc analysis indicated that VAF values significantly increased as the number of synergies increased from four to seven (p < 0.05). Also, with the same number of synergies, trial-specific synergy vector weights provided the highest VAF values, while subject-specific synergy vector weights provided the lowest VAF values for both S1 and S2 (Table 3). Overall, when using trial-specific synergy vector weights, measured muscle excitations were predicted with >95% VAF with four synergies for S1 and S2. When synergy vector weights were shared within the same speed, six synergies were needed for the left leg of S1, seven synergies were needed for the right leg of S1 and the left leg of S2, and four synergies were needed for the right leg of S2. As synergy vector weights were held constant across all trials, seven synergies were required for all legs except for the right leg of S2, where only 5 synergies were needed.

TABLE 3
www.frontiersin.org

TABLE 3. Mean ± standard deviation of VAF values for the reconstruction of measured muscle excitations with muscle synergy analysis across both subjects and all calibration trials with four to seven synergies when unmeasured synergy vector weights were assumed to be trial-specific, speed specific, or subject-specific.

3.2 First stage: Influence of methodological choices on SynX performance

In the first stage, the impact of different methodological combinations (i.e., number of synergies, category of unmeasured and residual synergy vector weights) on SynX performance was evaluated using r and RMSE values between predicted and measured iliopsoas muscle excitations. The three-factor ANOVA analyses revealed that the category of unmeasured synergy vector weights, category of residual synergy vector weights, and the number of synergies significantly affected both r and RMSE values that characterized SynX performance (p < 0.05). Several general observations were made by assessing the results for both psoas and iliacus and both subjects as a whole (Figure 2 and Supplementary Figure S2). First, without residual excitations being calibrated through optimization (labeled as “None” from “Params + SynX” calibration), r values were significantly lower, and RMSE values were significantly higher than those for any of the three categories of residual synergy vector weights (p < 0.05). Second, five and six synergies offered substantially higher r values and lower RMSE values than did four and seven synergies (p < 0.05), while no significant difference was detected between five and six synergies (p = 0.096). r values reached a maximum value, and RMSE values reached a minimum value at five synergies for the right leg of S2 or six synergies for both legs of S1 and the left leg of S2. Third, using speed-specific and trial-specific residual synergy vector weights yielded significantly better SynX performance than using subject-specific residual synergy vector weights (p < 0.05). Furthermore, when using five or six synergies, particularly trial-specific residual synergy vector weights, higher r and lower RMSE values were achieved when assuming them speed-specific. Fourth, with five or six synergies, using trial-specific synergy vector weights for the reconstruction of unmeasured muscle excitations provided the greatest r values and the lowest RMSE values, while subject-specific synergy vector weights provided the smallest r values and the biggest RMSE values (p < 0.05).

FIGURE 2
www.frontiersin.org

FIGURE 2. Synergy extrapolation performance for different methodological combinations when using the proposed EMG-driven calibration framework. (A) Pearson correlation coefficient r values and (B) root mean square error (RMSE) values for different methods of reconstructing psoas muscle activations across all calibration trials. (C) Pearson correlation coefficient r values and (D) root mean square error (RMSE) values for different methods of reconstructing iliacus muscle activations across all calibration trials. Unmeasured (bottom) and residual (side) synergy vector weights were categorized as either trial-specific, speed-specific, or subject-specific. Within the columns for 6 synergies, Orange boxes indicate the best SynX methodological combination (trial-specific unmeasured and speed-specific residual synergy vector weights) for calibration. The purple boxes indicate the best SynX methodological combination (speed-specific unmeasured and speed-specific residual synergy vector weights) for validation. Residual synergy vector weights categorized as “None” indicate calibration results for “Params + SynX” where no residual muscle excitations were predicted. These results suggest that “Params + SynX” calibration should be rejected due to unacceptable SynX performance.

Two “most reliable” SynX methodological combinations were identified for “calibration” and “validation” situations, respectively. As highlighted by oranges boxes in Figure 2 and Supplementary Figure S2, best SynX performance for both unmeasured muscled and both subjects was provided with six synergies, trial-specific unmeasured synergy vector weights, and speed-specific residual synergy vector weights, which could be used to analyze walking trials used in the calibration process. Further results obtained in the second stage for the “calibration” situation were described in Second stage: Evaluation of “calibration” situation of the results. As indicated by purple boxes in Figure 2 and Supplementary Figure S2, to analyze experimental walking trials not used in the calibration process or to generate computationally predicted walking motions, a SynX-performing methodological combination was chosen using six synergies and speed-specific synergy vector weights for the reconstruction of both unmeasured muscle excitations and residual muscle excitations. Further results obtained in the third stage for the “validation” situation were described in Third stage: Evaluation of “validation” situation of the results.

3.3 Second stage: Evaluation of “calibration” situation

For “Params + SynX + Res” calibration, the most reliable estimates of unmeasured muscle excitations and activations were achieved using six synergies, trial-specific unmeasured synergy vector weights, and speed-specific residual synergy vector weights (Figure 3). Compared to the results from the “Params” calibration, unmeasured muscle excitations (left two columns) were predicted with strong correlations for the left leg of S1 and both legs of S2 and with moderate correlations for the right leg of S1. The corresponding unmeasured muscle activations (right two columns) were strongly correlated with those from “Params” calibrations for both legs of S1 and the left leg of S2 and moderately correlated for the right leg of S2. SynX-estimated muscle excitations and activations were predicted with low RMSE values (<0.071) across all legs.

FIGURE 3
www.frontiersin.org

FIGURE 3. SynX-predicted muscle excitations and activations for psoas constructed with the best methodological combinations for the “calibration” situations (trial-specific unmeasured and speed-specific residual synergy vector weights with six synergies). Lines represent mean curves across calibration trials, and shaded areas represent ±1 standard deviation. r and RMSE values for muscle excitations and activations were calculated between “Params + SynX + Res” calibration (in pink) and “Params” calibration (in blue). “Params” calibration was performed using a complete set of EMG signals, where no muscle excitations were predicted. Data are reported for the complete gait cycle, where 0% indicates initial heel-strike and 100% indicates subsequent heel-strike for each leg of both subjects (right leg: paretic, left leg: nonparetic). SynX-predicted muscle excitations and activations for iliacus were of similar accuracy.

Compared with MAE values between model-predicted and experimental ID joint moments from “Params” calibrations, the MAE values were slightly lower but statistically comparable (p > 0.05) for all joints from “Params + SynX + Res” calibrations when residual excitations were calibrated but not used for joint moment calculation (Figure 4). However, the joint moment errors were significantly lowered (p < 0.05) when including residual muscle excitations in the prediction process (Supplementary Table S2 in the Appendix). On average, measured muscle activations generated from “Params + SynX + Res” calibration remained close to the those from “Params” calibration (Figure 5), where the greatest deviations were observed for the muscles that spanned the hip joint (e.g., glmed2). In addition, for measured and unmeasured muscles, muscle forces estimated from “Params + SynX + Res” calibration were in excellent agreement with those estimated from “Params” calibration in terms of both shape and magnitude (Figure 6). Moreover, in general, EMG-driven model parameter values were similar between “Params + SynX + Res” calibration and “Params” calibration (Figure 7). However, when additional variables (i.e., unmeasured muscle synergy vector weights and residual synergy vector weights) were tuned simultaneously, the pattern defined by the parameter magnitudes over all muscles was still retained for each model parameter.

FIGURE 4
www.frontiersin.org

FIGURE 4. Average joint moments across calibration trials and subjects from inverse dynamics (in red), “Params” calibration (in blue), and “Params + SynX + Res” calibration (in yellow). The results for “Params + SynX + Res” calibration were generated with the best methodological combination for calibration conditions (trial-specific unmeasured and speed-specific residual synergy vector weights with 6 synergies). For “Params + SynX + Res,” residual excitations were calibrated to improve the prediction of unmeasured muscle excitations but not used to calculate joint moments. Data are reported for the complete gait cycle, where 0% indicates initial heel-strike and 100% indicates subsequent heel-strike.

FIGURE 5
www.frontiersin.org

FIGURE 5. Average muscle activations across calibration trials and subjects from “Params” calibration (in blue), and “Params + SynX + Res” calibration (in yellow). Results for “Params + SynX + Res” calibration were generated using the best methodological combinations for analyzing experimentally measured walking motions (trial-specific unmeasured and speed-specific residual synergy vector weights with six synergies). Here, residual excitations were calibrated but not used to calculate muscle activations for “Params + SynX + Res.” Data are reported for the complete gait cycle, where 0% indicates initial heel-strike and 100% indicates subsequent heel-strike.

FIGURE 6
www.frontiersin.org

FIGURE 6. Average muscle forces across calibration trials and subjects from “Params” calibration (in blue) and “Params + SynX + Res” calibration (in yellow). Results for “Params + SynX + Res” calibration were produced using the best methodological combinations for analyzing experimentally measured walking motions (trial-specific unmeasured and speed-specific residual synergy vector weights with six synergies). Here, residual excitations were calibrated but not used to calculate muscle forces for “Params + SynX + Res.” Data are reported for the complete gait cycle, where 0% indicates initial heel-strike and 100% indicates subsequent heel-strike.

FIGURE 7
www.frontiersin.org

FIGURE 7. Average EMG-driven model parameters of both subjects from “Params” calibration (in blue) and “Params + SynX + Res” calibration (in orange) Results for “Params + SynX + Res” calibration were produced using the best performing methodological combinations for analyzing experimentally measured walking motions (trial-specific unmeasured and speed-specific residual synergy vector weights with six synergies). The upper and lower bounds for each of the four activation dynamics model parameters during optimization have been indicated by green vertical dash-sot lines, where the upper and lower bounds for the scaling factors of optimal fiber lengths and tendon slack lengths were [0.6, 1.4] for all muscles.

3.4 Third stage: Evaluation of “validation” situation

Best SynX methodological choices to analyze experimental walking trials not used in the calibration process or to generate computationally predicted walking motions was a combination of speed-specific unmeasured synergy vector weights speed-specific residual vector weights and with 6 synergies (Figure 2). From “Params + SynX + Res” calibrations with this method, estimated muscle excitations for “iliopsoas” (left two columns) were correlated with the results from “Params” calibration strongly for left legs of both subjects and moderately for right legs of both subjects (Figure 8). There was strong correlations between muscle activations for “psoas” (right two columns) from “Params + SynX + Res” calibrations and the ones from “Params” calibrations for both legs of S1 and left leg of S2, while the correlation was moderate for right leg of S2. As for magnitude, the highest errors (RMSE = 0.1) occurred in the muscle activations for right leg of S2.

FIGURE 8
www.frontiersin.org

FIGURE 8. SynX-predicted muscle excitations and activations for psoas constructed with the best methodological combinations for the “validation” situations (speed-specific unmeasured and speed-specific residual synergy vector weights with six synergies). Lines represent mean curves across calibration trials, and shaded areas represent ±1 standard deviation. r and RMSE values for muscle excitations and activations were calculated between “Params + SynX + Res” calibration (in pink) and “Params” calibration (in blue). “Params” calibration was performed using a complete set of EMG signals, where no muscle excitations were predicted. Data are reported for the complete gait cycle, where 0% indicates initial heel-strike and 100% indicates subsequent heel-strike for each leg of both subjects (right leg: paretic, left leg: nonparetic). SynX-predicted muscle excitations and activations for iliacus were of similar accuracy (not shown).

“Params + SynX + Res” calibrations with this methodological combination were able to capture experimental joint moments with reasonable accuracy (Table 4 and Figure 9A). With respect to the range of variation assumed by the average experimental joint moments, average MAE values were 0.09 for HipFE, 0.10 for HipAA, 0.22 for HipRot, 0.11 for KneeFE, 0.07 for AnklePD and 0.21 for AnkleIE. Moreover, calibrated EMG-driven models and speed-specific unmeasured synergy vector weights were capable of predicting prediction joint moments for validation trials with comparable accuracy (Table 4 and Figure 9B). Compared with the joint moments estimated from “Params”, “Params + SynX + Res” calibrations matched the ID joint moments slightly better at three hip DOFs and slightly worse at knee DOF and two ankle DOFs for both calibration trials and validation trials, whereas the difference were not statistically significant (p = 0.092, Table 4 and Figure 9), Additionally, for both calibration and validation trials, “Params + SynX + Res” calibration using this methodological combination could provide similar estimates of muscle activations (Supplementary Figure S5) or muscle forces (Supplementary Figure S6) to the results predicted from “Params” calibration case in terms of both shape and magnitude.

TABLE 4
www.frontiersin.org

TABLE 4. Mean absolute error (MAE) values calculated between joint moments found from inverse dynamics and either “Params” calibration or “Params + SynX + Res” calibration.

FIGURE 9
www.frontiersin.org

FIGURE 9. Average joint moments across trials and subjects for (A) calibration trials and (B) validation trials from inverse dynamics (in red), “Params” calibration (in blue), and “Params + SynX + Res” calibration (in yellow). Results for “Params + SynX + Res” calibration were produced using the best methodological combinations for generating computationally predicted walking motions (speed-specific unmeasured and residual synergy vector weights with six synergies). Here, residual excitations were calibrated but not used to calculate joint moments for “Params + SynX + Res.” Data are reported for the complete gait cycle, where 0% indicates initial heel-strike and 100% indicates subsequent heel-strike.

4 Discussion

In this paper, a novel computational approach was presented to address the problem of a small number of missing EMG signals during EMG-driven model calibration (Figure 1). SynX was employed to predict unmeasured muscle excitations by linearly combining synergy excitations extracted from measured muscle excitations, where the contribution of each synergy excitations to an unmeasured muscle excitation was determined by a set of unmeasured synergy vector weights. Meanwhile, residual muscle excitations constructed as a linear combination of measured synergy excitations were added to measured muscle excitations. We explored the influence of several important SynX methodological choices on prediction accuracy of unmeasured muscle excitations and muscle activations. We then identified two SynX methodological combinations to analyze experimental walking trials used in the calibration process and not used in the calibration process. First, a combination of trial-specific unmeasured synergy vector weights, speed-specific residual synergy vector weights and six synergies should be hired for analyzing experimental walking motions included in the calibration process. This methodological combination can not only provide accurate and reliable missing muscle excitations and activations, but estimates of muscle forces and joint moments that stayed close to the ones obtained from “Params” calibration for both subjects. Second, with a combination of speed-specific unmeasured synergy vector weights, speed-specific residual synergy vector weights and six synergies, we were still capable of estimating these important biomechanical quantities with reasonable accuracy. More importantly, these methods permitted us to analyze experimental walking trials not used in the calibration process or generate computationally predicted walking motions.

Synergy extrapolation (SynX) has been an emerging approach to predict unmeasured muscle excitations by the use of muscle synergy information extracted from measured muscle excitations (Bianco et al., 2018; Ao et al., 2020). Bianco et al. investigated the theoretical feasibility of using synergy excitations extracted from a group of eight “included” muscle excitations treated as measured to construct muscle excitations for a group of eight “excluded” muscle excitations treated as unmeasured (Bianco et al., 2018). To take a step forward, a previous paper from the author has shown that SynX could provide accurate estimates of unmeasured muscle excitations, where unmeasured synergy vector weights were identified through non-linear optimizations by incorporating a well-calibrated EMG-driven model (Ao et al., 2020). However, the “well-calibrated” EMG-driven model was obtained with the knowledge of muscle excitations being treated as unmeasured. Thus, this paper had achieved to implement SynX without demands of EMG data that were sought to be predict, which marked the feasibility of SynX in practical applications. Moreover, EMG-driven model calibration with SynX was formulated such that both unmeasured muscle excitations and EMG-driven model parameters were defined simultaneously, which provided advantages for both EMG-driven model personalization and missing muscle excitation prediction.

Calibration of synergy-structured residual muscle excitations were included in the proposed framework to improve the accuracy of SynX-predicted unmeasured muscle excitations (Figure 1). Noted that even though we performed a sequence of steps to improve model correctness and subject-specificity prior to EMG-driven model calibration, there were still some resources of errors that affected accuracy of joint moment estimation, such as marker location errors, soft tissue movement errors, and surface EMG measurement errors (Lloyd and Besier, 2003; Manal and Buchanan, 2013; Meyer et al., 2017). Therefore, when unmeasured muscle excitations became design variables that were iteratively adjusted within EMG-driven model calibration for SynX, they inclined to deviate from experimental muscle excitations such that joint moment matching errors were minimized (see “Params + SynX” calibrated results in Figure 2 and Supplementary Figure S2 which were labeled as “None” for the category of residual synergy vector weights). Residual muscle excitations were introduced to account for joint moment matching errors, which facilitated to prevent predicted missing muscle excitations from excessively compensating for joint moment prediction inaccuracy through optimization and becoming inaccurate as a consequence.

Even though residual muscle excitations were needed to predict within “Params + SynX + Res” optimization, our EMG-driven model calibration method could still generate two sets of joint moments, depending on whether residual muscle excitations were included or not for joint moment calculation. To illustrate the discrepancy in results between “Params + SynX + Res” and “Params” calibrations, we consistently presented the muscle activations (Figure 5), muscle forces (Figure 6) and joint moments (Figure 4) calculated without inclusion of residual muscle excitations. However, joint moment prediction errors could be substantially reduced (by averagely 3 Nm) when including residual muscle excitations for joint moment calculation (Supplementary Table S2). Thus, joint moments outputted from our framework could be highly dependent on the actual demands, where residual muscle excitations could be included in the prediction of joint moments, if the difference between model-predicted and experimental joint moments are required to be as low as possible.

As a follow-up step of “Params” calibration, we calibrated the residual excitations needed for all the muscles (including psoas and iliacus) to match joint moments from inverse dynamics better (termed “Res” calibration), where the EMG-driven model parameters were held constant with the values obtained through “Params” calibration. We found that joint moment tracking accuracy was significantly improved by including synergy-structured residual excitations (p < 0.05), where the extent of improvements were dependent on the category of residual synergy vector weights (Supplementary Figure S3). More interestingly, we observed that within the regions of gait cycle where SynX-predicted muscle excitations and activations from “Params + SynX + Res” calibration could not match well with the ones from “Params” calibration, they showed good agreement with the ones from “Res” calibrations (Supplementary Figures S8, S9). Therefore, we believed that SynX-predicted muscle excitations and activations within our framework were likely superposition of “experimental” muscle activities and corresponding synergy-structured residual muscle activities needed for psoas and iliacus to improve joint moment matching.

For several important reasons, synergy structures were imposed on both unmeasured muscle excitations and residual muscle excitations in our framework. First, unlike SO that solved a time frame of muscle activation at a time (Anderson and Pandy, 2001; Ackermann and van den Bogert, 2010; Sartori et al., 2014; Zonnino and Sergi, 2019), unmeasured and residual synergy vector weights were time-invariant, which opened the possibility of a single-layer optimization process to simultaneously achieve EMG-driven model personalization and missing muscle excitation prediction. Second, muscle excitation-activation relationships were defined by a set of differential equations (Meyer et al., 2017), synergy-structured excitations could be resolved over all time frames at a time, which allowed us to perform reconstruction at level of muscle excitations and thus to decrease the search space for the optimization in comparison with SO-based approaches. Third, the problem of finding unknown time-varying activations was reduced to a problem of finding a smaller number of unmeasured and residual synergy vector weights, which significantly decreased the search space for the optimization in comparison with SO-based approaches. Fourth, due to inherent constraints of dependence between time frames, synergy-structured residual excitations showed predictably worse joint moment matching performance than did SO-based residual activations (Supplementary Figure S3). However, they could still lower joint moment matching errors by significant amount, especially when using trial-specific synergy vector weights (Supplementary Figure S3). When both unmeasured muscle excitations and residual excitations were consistently constructed using synergy concepts, the joint moment matching errors that resultant unmeasured muscle excitations could potentially over-compensate for would be sufficiently accounted for by the resultant residual excitations.

An optimization problem was formulated within our EMG-driven modeling method that minimized four primary cost terms (Eq. 15) simultaneously, where a combination of maximum allowable deviation (MAD) values needed for normalization of cost values was determined through a series of sensitivity tests such that satisfactory SynX outcomes were consistently generated across legs of both subjects (Supplementary Figure S1). On the one hand, joint moment tracking errors were minimized when residual muscle excitations were (Eq. 16) and were not (Eq. 17) used for joint moment calculation. Both terms were included such that residual excitations for measured muscles were estimated with reasonable accuracy and EMG-driven model parameters for measured muscles stayed as close as possible to the ones derived from “Params” calibration. On the other hand, unmeasured and residual muscle activations were minimized such that the optimization did not incline to converge to a set of overestimated solutions. Here, we implemented minimization at the activation-level instead of excitation-level, because muscle activations were more consistent in terms of magnitudes across muscles and subjects, which facilitated finding a consistent set of maximum allowable deviation values (MAD) that could provide acceptable outcomes across legs of both subjects.

There were a number of methodological choices required to make during muscle synergy analysis (MSA) that could influence the results of extracted measured synergy excitations and corresponding synergy extrapolation performance, such as EMG normalization methods, number of synergies, matrix decomposition algorithm (Ivanenko et al., 2005; Tresch et al., 2006; Steele et al., 2013; Oliveira et al., 2014; Banks et al., 2017; Shuman et al., 2017; Ebied et al., 2018; Gallina et al., 2018; Ao et al., 2020) In a previous study (Ao et al., 2020), the author reported that PCA provided more accurate, reliable, and efficient estimates of unmeasured muscle excitations than NMF when using a well-calibrated EMG-driven model. A similar observation was made in this study that when the number of synergies and categories of unmeasured and residual synergy vector weights were matched, NMF could not produce estimates of unmeasured muscle excitations and activations from “Params” calibration as accurately as was PCA (Supplementary Figure S7). As highlighted in the previous paper, the reasons PCA outperformed NMF for SynX could be due to nonnegativity constraints for NMF and extra design variables for PCA, both of which could make the feasible search space of NMF more restricted than that of PCA. Beyond these, PCA was especially advantageous in our framework because it allowed residual excitations to be both positive and negative, which could be beneficial to achievement of lower joint moment errors.

The principal methodological decisions explored in this study were the number of synergies and variability of synergy vector weights across trials. With an increasing number of synergies, SynX performance exhibited non-monotonic behavior for both subjects, where six synergies provided the generally best SynX performance and EMG-driven model calibration outcomes for both legs of S1 and left leg of S2 and five synergies provided the best overall results for right leg of S2. The results were in great agreement with the findings reported in (Ao et al., 2020). Additionally, based on the assumptions about the variability of synergy vector weight across walking trials, we categorized them associated with unmeasured and residual muscle excitations as trial-specific, speed-specific, and subject-specific, respectively, while different concatenation strategies were used to extract corresponding synergy excitations. We evaluated the results with all possible methodological combinations (Figure 2 and Supplementary Figure S2), and it was indicated that with matched number of synergies, trial-specific unmeasured synergy vector weights and speed-specific residual synergy vector weights generated the best SynX performance for most legs. For both measured (Table 3) and unmeasured (Figure 2 and Supplementary Figure S2) muscle excitations, a small number of synergies (e.g., <5) or subject-specific synergy vector weights may not be sufficient to account for their variance. As the number of synergies or the variability of synergy vector weights across trials increased, the additional degrees of freedom in the optimization allowed the optimizer to lower the joint moment matching errors (Supplementary Figure S4). However, when the flexibility determined by the number of synergies and category of synergy vector weights was above a certain level, the joint moment tracking errors dropped below those achieved by calibration with a complete set of EMG data (“Params” calibration), where the prediction of unmeasured muscle excitations and activations becoming less accurate as a consequence. Another important merit of using speed-specific residual synergy vector weights was that computational cost would be drastically reduced because we normally have way more channels of measured EMGs than channels of unmeasured EMGs when using our framework.

Even though we suggested employing six synergies, which provided the most accurate predictions for most legs, the overall SynX and model calibration performance peaked at five synergies for the right leg of S2 (Figure 2 and Supplementary Figure S2). The underlying reason for this inconsistency could be that subject S2 was a low-functioning post-stroke subject, whose right leg was paretic. Clark et al. reported that post-stroke subjects had fewer muscle synergies which could result from the merging of the synergies observed in healthy controls due to impaired locomotor coordination and the reduced independence of neural control signals (Clark et al., 2009). We made a consistent observation within our MSA results that the right leg of S2 had generally higher VAF values when the number of synergies and category of synergy vector weights were matched. Also, when assuming synergy vector weights speed-specific and subject-specific, the right leg of S2 required four or five synergies to account for over 95% of the variability in measured muscle excitations, while the left leg of S2 and both legs of S1 needed six or seven synergies. Therefore, it is reasonable to speculate that for the paretic leg of S2, five synergies were sufficient for SynX when hiring our framework, while six synergies could introduce additional degrees of freedom to the optimization problem and influence the SynX performance adversely.

This study possessed several limitations that could help inform future research efforts. First, this study validated the performance of our framework using gait datasets from two subjects post-stroke who had a complete set of EMG data from muscles in the lower extremities. Broader investigations need to be done in diverse subject populations of larger sample sizes. Second, our framework was validated in cases where EMG data from only one important hip muscle group (iliopsoas) was assumed to be missing at a time. The approach should be evaluated more extensively in cases where more important deep muscles are assumed to be unmeasured individually or concurrently. Third, we developed the framework using walking data with two representative speeds. We are planning to explore the framework’s feasibility for other types of dynamic movement conditions and experimental scenarios, such as stair climbing and running. Lastly, to achieve the purpose of computational prediction for the trials withheld from calibration, we suggested assuming unmeasured and residual synergy vector weights speed-specific for walking. When our framework is applied for diverse motion tasks, it is worthwhile exploring the feasibility of making a hybrid assumption about the variability of synergy vector weights, such as speed-specific synergy vector weights mixed with task-specific synergy vector weights.

5 Conclusion

In conclusion, we have presented a novel computational method that can address the problem of a small number of missing EMG signals during EMG-driven model calibration. Within this method, synergy excitations extracted from measured muscle excitations using PCA were linearly combined to reconstruct unmeasured muscle excitations (termed “synergy extrapolation” or “SynX”) and residual muscle excitations for measured muscles. The unmeasured and residual synergy vector weights defining the contribution of each measured synergy excitation to all muscle excitations were identified simultaneously with EMG-driven model parameters through a multi-objective optimization. The study also assessed how SynX methodological choices (i.e., number of synergies and category of unmeasured synergy vector weights and residual synergy vector weights) influenced SynX performance. By comparing with results from EMG-driven model calibration using a complete set of EMG signals, we identified two SynX methodological combinations for the purposes of analyzing experimental walking trials used in the calibration process and not used in the calibration process, respectively. Both methodological combinations consistently provided accurate, reliable, and efficient estimates of both SynX-relevant quantities (i.e., missing muscle excitations and activations) and biomechanical variables (i.e., muscle forces and joint moments). This computational approach opens up possibilities for the personalization of EMG-driven musculoskeletal models when difficulties exist with collecting EMG signals from important muscles.

Data availability statement

The experimental data, OpenSim model, and Matlab code used to perform this study are available at https://simtk.org/projects/synx.

Ethics statement

All experimental procedures were approved by the University of Florida Health Science Center Institutional Review Board (IRB-01), and the subject provided written informed consent prior to participation.

Author contributions

CP and BF designed and performed the experiments; DA and MV wrote the programs; DA analyzed the data, prepared figures, and drafted the manuscript; DA wrote the manuscript; MS, MV and BF revised the manuscript; DA, MV, MS, CP, and BF approved the final version of the manuscript.

Funding

This work was conducted with support from the Cancer Prevention and Research Institute of Texas (CPRIT) RR170026.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fbioe.2022.962959/full#supplementary-material

References

Ackermann, M., and van den Bogert, A. J. (2010). Optimality principles for model-based prediction of human gait. J. Biomech. 43, 1055–1060. doi:10.1016/j.jbiomech.2009.12.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Amarantini, D., and Martin, L. (2004). A method to combine numerical optimization and EMG data for the estimation of joint moments under dynamic conditions. J. Biomech. 37, 1393–1404. doi:10.1016/j.jbiomech.2003.12.020

PubMed Abstract | CrossRef Full Text | Google Scholar

An, K. N., Takahashi, K., Harrigan, T. P., and Chao, E. Y. (1984). Determination of muscle orientations and moment arms. J. Biomech. Eng. 106, 280–282. doi:10.1115/1.3138494

PubMed Abstract | CrossRef Full Text | Google Scholar

Anderson, F. C., and Pandy, M. G. (2001). Static and dynamic optimization solutions for gait are practically equivalent. J. Biomech. 34, 153–161. doi:10.1016/s0021-9290(00)00155-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Ao, D., Shourijeh, M. S., Patten, C., and Fregly, B. J. (2020). Evaluation of synergy extrapolation for predicting unmeasured muscle excitations from measured muscle synergies. Front. Comput. Neurosci. 108, 588943. doi:10.3389/fncom.2020.588943

PubMed Abstract | CrossRef Full Text | Google Scholar

Banks, C. L., Pai, M. M., McGuirk, T. E., Fregly, B. J., and Patten, C. (2017). Methodological choices in muscle synergy analysis impact differentiation of physiological characteristics following stroke. Front. Comput. Neurosci. 11, 78. doi:10.3389/fncom.2017.00078

PubMed Abstract | CrossRef Full Text | Google Scholar

Bianco, N. A., Patten, C., and Fregly, B. J. (2018). Can measured synergy excitations accurately construct unmeasured muscle excitations? J. Biomech. Eng. 140, 11011. doi:10.1115/1.4038199

CrossRef Full Text | Google Scholar

Buchanan, T. S., Lloyd, D. G., Manal, K., and Besier, T. F. (2005). Estimation of muscle forces and joint moments using a forward-inverse dynamics model. Med. Sci. Sports Exerc. 37, 1911–1916. doi:10.1249/01.mss.0000176684.24008.6f

PubMed Abstract | CrossRef Full Text | Google Scholar

Cano-de-la-Cuerda, R., Pérez-de-Heredia, M., Miangolarra-Page, J. C., Munoz-Hellín, E., and Fernández-de-Las-Peñas, C. (2010). Is there muscular weakness in Parkinson’s disease? Am. J. Phys. Med. Rehabil. 89, 70–76. doi:10.1097/phm.0b013e3181a9ed9b

PubMed Abstract | CrossRef Full Text | Google Scholar

Cappellini, G., Ivanenko, Y. P., Poppele, R. E., and Lacquaniti, F. (2006). Motor patterns in human walking and running. J. Neurophysiol. 95, 3426–3437. doi:10.1152/jn.00081.2006

PubMed Abstract | CrossRef Full Text | Google Scholar

Clark, D. J., Ting, L. H., Zajac, F. E., Neptune, R. R., and Kautz, S. A. (2009). Merging of healthy motor modules predicts reduced locomotor performance and muscle coordination complexity post-stroke. J. Neurophysiol. 103, 844–857. doi:10.1152/jn.00825.2009

PubMed Abstract | CrossRef Full Text | Google Scholar

Contessa, P., and Luca, C. J. (2013). Neural control of muscle force: Indications from a simulation model. J. Neurophysiol. 109, 1548–1570. doi:10.1152/jn.00237.2012

PubMed Abstract | CrossRef Full Text | Google Scholar

Correa, T. A., Crossley, K. M., Kim, H. J., and Pandy, M. G. (2010). Contributions of individual muscles to hip joint contact force in normal walking. J. Biomech. 43, 1618–1622. doi:10.1016/j.jbiomech.2010.02.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Crowninshield, R. D., and Brand, R. A. (1981). A physiologically based criterion of muscle force prediction in locomotion. J. Biomech. 14, 793–801. doi:10.1016/0021-9290(81)90035-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Del Vecchio, A., Ubeda, A., Sartori, M., Azorin, J. M., Felici, F., and Farina, D. (2018). Central nervous system modulates the neuromechanical delay in a broad range for the control of muscle force. J. Appl. Physiol. (1985). 125, 1404–1410. doi:10.1152/japplphysiol.00135.2018

PubMed Abstract | CrossRef Full Text | Google Scholar

Delp, S. L., Anderson, F. C., Arnold, A. S., Loan, P., Habib, A., John, C. T., et al. (2007). OpenSim: Open-source software to create and analyze dynamic simulations of movement. IEEE Trans. Biomed. Eng. 54, 1940–1950. doi:10.1109/tbme.2007.901024

PubMed Abstract | CrossRef Full Text | Google Scholar

Ebied, A., Kinney-Lang, E., Spyrou, L., and Escudero, J. (2018). Evaluation of matrix factorisation approaches for muscle synergy extraction. Med. Eng. Phys. 57, 51–60. doi:10.1016/j.medengphy.2018.04.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Gallina, A., Garland, S. J., and Wakeling, J. M. (2018). Identification of regional activation by factorization of high-density surface EMG signals: A comparison of principal component analysis and non-negative matrix factorization. J. Electromyogr. Kinesiol. 41, 116–123. doi:10.1016/j.jelekin.2018.05.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Gurchiek, R. D., Ursiny, A. T., and McGinnis, R. S. (2020). A Gaussian process model of muscle synergy functions for estimating unmeasured muscle excitations using a measured subset. IEEE Trans. Neural Syst. Rehabil. Eng. 28, 2478–2487. doi:10.1109/tnsre.2020.3028052

PubMed Abstract | CrossRef Full Text | Google Scholar

He, J., Levine, W. S., and Loeb, G. E. (1991). Feedback gains for correcting small perturbations to standing posture. IEEE Trans. Autom. Contr. 36, 322–332. doi:10.1109/9.73565

CrossRef Full Text | Google Scholar

Heintz, S., and Gutierrez-Farewik, E. M. (2007). Static optimization of muscle forces during gait in comparison to EMG-to-force processing approach. Gait Posture 26, 279–288. doi:10.1016/j.gaitpost.2006.09.074

PubMed Abstract | CrossRef Full Text | Google Scholar

Herzog, W., and Binding, P. (1992). Predictions of antagonistic muscular activity using nonlinear optimization. Math. Biosci. 111, 217–229. doi:10.1016/0025-5564(92)90071-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Hill, A. V. (1938). The heat of shortening and the dynamic constants of muscle. Proc. R. Soc. Lond. Ser. B-Biological Sci. 126, 136–195.

Google Scholar

Hoang, H. X., Diamond, L. E., Lloyd, D. G., and Pizzolato, C. (2019). A calibrated EMG-informed neuromusculoskeletal model can appropriately account for muscle co-contraction in the estimation of hip joint contact forces in people with hip osteoarthritis. J. Biomech. 83, 134–142. doi:10.1016/j.jbiomech.2018.11.042

PubMed Abstract | CrossRef Full Text | Google Scholar

Ivanenko, Y. P., Cappellini, G., Dominici, N., Poppele, R. E., and Lacquaniti, F. (2005). Coordination of locomotion with voluntary movements in humans. J. Neurosci. 25, 7238–7253. doi:10.1523/jneurosci.1327-05.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, H. J., Fernandez, J. W., Akbarshahi, M., Walter, J. P., Fregly, B. J., and Pandy, M. G. (2009). Evaluation of predicted knee‐joint muscle forces during gait using an instrumented knee implant. J. Orthop. Res. 27, 1326–1331. doi:10.1002/jor.20876

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, G., Shourijeh, M. S., Ao, D., Patten, C., and Fregly, B. J. (2020). How well do commonly used Co-contraction indices approximate lower limb joint stiffness trends during gait? bioRxiv 8, 588908. doi:10.3389/fbioe.2020.588908

PubMed Abstract | CrossRef Full Text | Google Scholar

Lloyd, D. G., and Besier, T. F. (2003). An EMG-driven musculoskeletal model to estimate muscle forces and knee joint moments in vivo. J. Biomech. 36, 765–776. doi:10.1016/s0021-9290(03)00010-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Manal, K., and Buchanan, T. S. (2003). A one-parameter neural activation to muscle activation model: Estimating isometric joint moments from electromyograms. J. Biomech. 36, 1197–1202. doi:10.1016/s0021-9290(03)00152-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Manal, K., and Buchanan, T. S. (2013). An electromyogram-driven musculoskeletal model of the knee to predict in vivo joint contact forces during normal and novel gait patterns. J. Biomech. Eng. 135, 021014. doi:10.1115/1.4023457

PubMed Abstract | CrossRef Full Text | Google Scholar

Martin, J. A., Brandon, S. C. E., Keuler, E. M., Hermus, J. R., Ehlers, A. C., Segalman, D. J., et al. (2018). Gauging force by tapping tendons. Nat. Commun. 9, 1592–1599. doi:10.1038/s41467-018-03797-6

PubMed Abstract | CrossRef Full Text | Google Scholar

McLean, S. G., Huang, X., and van den Bogert, A. J. (2005). Association between lower extremity posture at contact and peak knee valgus moment during sidestepping: Implications for ACL injury. Clin. Biomech. (Bristol, Avon. 20, 863–870. doi:10.1016/j.clinbiomech.2005.05.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Menegaldo, L. L., de Oliveira, L. F., and Minato, K. K. (2014). EMGD-FE: An open source graphical user interface for estimating isometric muscle forces in the lower limb using an EMG-driven model. Biomed. Eng. Online 13, 37–12. doi:10.1186/1475-925x-13-37

PubMed Abstract | CrossRef Full Text | Google Scholar

Menegaldo, L. L., de Toledo Fleury, A., and Weber, H. I. (2004). Moment arms and musculotendon lengths estimation for a three-dimensional lower-limb model. J. Biomech. 37, 1447–1453. doi:10.1016/j.jbiomech.2003.12.017

PubMed Abstract | CrossRef Full Text | Google Scholar

Meyer, A. J., Eskinazi, I., Jackson, J. N., Rao, A. V., Patten, C., and Fregly, B. J. (2016). Muscle synergies facilitate computational prediction of subject-specific walking motions. Front. Bioeng. Biotechnol. 4, 77. doi:10.3389/fbioe.2016.00077

PubMed Abstract | CrossRef Full Text | Google Scholar

Meyer, A. J., Patten, C., and Fregly, B. J. (2017). Lower extremity EMG-driven modeling of walking with automated adjustment of musculoskeletal geometry. PLoS One 12, e0179698. doi:10.1371/journal.pone.0179698

PubMed Abstract | CrossRef Full Text | Google Scholar

Michaud, F., Shourijeh, M. S., Fregly, B. J., and Cuadrado, J. (2020). Do muscle synergies improve optimization prediction of muscle activations during gait? Front. Comput. Neurosci. 14, 54. doi:10.3389/fncom.2020.00054

PubMed Abstract | CrossRef Full Text | Google Scholar

Neptune, R. R., Clark, D. J., and Kautz, S. A. (2009). Modular control of human walking: A simulation study. J. Biomech. 42, 1282–1287. doi:10.1016/j.jbiomech.2009.03.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Oliveira, A. S., Gizzi, L., Farina, D., and Kersting, U. G. (2014). Motor modules of human locomotion: Influence of EMG averaging, concatenation, and number of step cycles. Front. Hum. Neurosci. 8, 335. doi:10.3389/fnhum.2014.00335

PubMed Abstract | CrossRef Full Text | Google Scholar

Péter, A., Andersson, E., Hegyi, A., Finni, T., Tarassova, O., Cronin, N., et al. (2019). Comparing surface and fine-wire electromyography activity of lower leg muscles at different walking speeds. Front. Physiol. 10, 1283. doi:10.3389/fphys.2019.01283

PubMed Abstract | CrossRef Full Text | Google Scholar

Pizzolato, C., Lloyd, D. G., Sartori, M., Ceseracciu, E., Besier, T. F., Fregly, B. J., et al. (2015). Ceinms: A toolbox to investigate the influence of different neural control solutions on the prediction of muscle excitation and joint moments during dynamic motor tasks. J. Biomech. 48, 3929–3936. doi:10.1016/j.jbiomech.2015.09.021

PubMed Abstract | CrossRef Full Text | Google Scholar

Rajagopal, A., Dembia, C. L., DeMers, M. S., Delp, D. D., Hicks, J. L., and Delp, S. L. (2016). Full-body musculoskeletal model for muscle-driven simulation of human gait. IEEE Trans. Biomed. Eng. 63, 2068–2079. doi:10.1109/tbme.2016.2586891

PubMed Abstract | CrossRef Full Text | Google Scholar

Reinbolt, J. A., Schutte, J. F., Fregly, B. J., Koh, B. Il, Haftka, R. T., George, A. D., et al. (2005). Determination of patient-specific multi-joint kinematic models through two-level optimization. J. Biomech. 38, 621–626. doi:10.1016/j.jbiomech.2004.03.031

PubMed Abstract | CrossRef Full Text | Google Scholar

Sartori, M., Farina, D., and Lloyd, D. G. (2014). Hybrid neuromusculoskeletal modeling to best track joint moments using a balance between muscle excitations derived from electromyograms and optimization. J. Biomech. 47, 3613–3621. doi:10.1016/j.jbiomech.2014.10.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Sartori, M., Gizzi, L., Lloyd, D. G., and Farina, D. (2013). A musculoskeletal model of human locomotion driven by a low dimensional set of impulsive excitation primitives. Front. Comput. Neurosci. 7, 79. doi:10.3389/fncom.2013.00079

PubMed Abstract | CrossRef Full Text | Google Scholar

Sartori, M., Reggiani, M., Farina, D., and Lloyd, D. G. (2012). EMG-driven forward-dynamic estimation of muscle force and joint moment about multiple degrees of freedom in the human lower extremity. PLoS One 7, e52618. doi:10.1371/journal.pone.0052618

PubMed Abstract | CrossRef Full Text | Google Scholar

Sasaki, K., and Neptune, R. R. (2010). Individual muscle contributions to the axial knee joint contact force during normal walking. J. Biomech. 43, 2780–2784. doi:10.1016/j.jbiomech.2010.06.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Sauder, N. R., Meyer, A. J., Allen, J. L., Ting, L. H., Kesar, T. M., and Fregly, B. J. (2019). Computational design of FastFES treatment to improve propulsive force symmetry during post-stroke gait: A feasibility study. Front. Neurorobot. 13, 80. doi:10.3389/fnbot.2019.00080

PubMed Abstract | CrossRef Full Text | Google Scholar

Schellenberg, F., Oberhofer, K., Taylor, W. R., and Lorenzetti, S. (2015). Review of modelling techniques for in vivo muscle force estimation in the lower extremities during strength training. Comput. Math. Methods Med. 2015, 483921. doi:10.1155/2015/483921

PubMed Abstract | CrossRef Full Text | Google Scholar

Serrancolí, G., Kinney, A. L., Fregly, B. J., and Font-Llagunes, J. M. (2016). Neuromusculoskeletal model calibration significantly affects predicted knee contact forces for walking. J. Biomech. Eng. 138. doi:10.1115/1.4033673

PubMed Abstract | CrossRef Full Text | Google Scholar

Seth, A., Hicks, J. L., Uchida, T. K., Habib, A., Dembia, C. L., Dunne, J. J., et al. (2018). OpenSim: Simulating musculoskeletal dynamics and neuromuscular control to study human and animal movement. PLoS Comput. Biol. 14, e1006223. doi:10.1371/journal.pcbi.1006223

PubMed Abstract | CrossRef Full Text | Google Scholar

Shao, Q., Bassett, D. N., Manal, K., and Buchanan, T. S. (2009). An EMG-driven model to estimate muscle forces and joint moments in stroke patients. Comput. Biol. Med. 39, 1083–1088. doi:10.1016/j.compbiomed.2009.09.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Shourijeh, M. S., and Fregly, B. J. (2020). Muscle synergies modify optimization estimates of joint stiffness during walking. J. Biomech. Eng. 142, 011011. doi:10.1115/1.4044310

PubMed Abstract | CrossRef Full Text | Google Scholar

Shuman, B. R., Schwartz, M. H., and Steele, K. M. (2017). Electromyography data processing impacts muscle synergies during gait for unimpaired children and children with cerebral palsy. Front. Comput. Neurosci. 11, 50. doi:10.3389/fncom.2017.00050

PubMed Abstract | CrossRef Full Text | Google Scholar

Steele, K. M., Tresch, M. C., and Perreault, E. J. (2013). The number and choice of muscles impact the results of muscle synergy analyses. Front. Comput. Neurosci. 7, 105. doi:10.3389/fncom.2013.00105

PubMed Abstract | CrossRef Full Text | Google Scholar

Taylor, R. (1990). Interpretation of the correlation coefficient: A basic review. J. Diagn. Med. Sonogr. 6, 35–39. doi:10.1177/875647939000600106

CrossRef Full Text | Google Scholar

Ting, L. H., and Chvatal, S. A. (2010). Mot. Control theor. Exp. Appl. Oxf. New York: Univ. Press, 102v–138.Decomposing muscle activity in motor tasks

Google Scholar

Tresch, M. C., Cheung, V. C. K., and d’Avella, A. (2006). Matrix factorization algorithms for the identification of muscle synergies: Evaluation on simulated and experimental data sets. J. Neurophysiol. 95, 2199–2212. doi:10.1152/jn.00222.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

Vilimek, M. (2007). Musculotendon forces derived by different muscle models. Acta Bioeng. Biomech. 9, 41–47.

PubMed Abstract | Google Scholar

Walter, J. P., Kinney, A. L., Banks, S. A., D’Lima, D. D., Besier, T. F., Lloyd, D. G., et al. (2014). Muscle synergies may improve optimization prediction of knee contact forces during walking. J. Biomech. Eng. 136, 021031. doi:10.1115/1.4026428

PubMed Abstract | CrossRef Full Text | Google Scholar

Zajac, F. E. (1989). Muscle and tendon: Properties, models, scaling, and application to biomechanics and motor control. Crit. Rev. Biomed. Eng. 17, 359–411.

PubMed Abstract | Google Scholar

Zonnino, A., and Sergi, F. (2019). Model-based estimation of individual muscle force based on measurements of muscle activity in forearm muscles during isometric tasks. IEEE Trans. Biomed. Eng. 1, 134–145. doi:10.1109/tbme.2019.2909171

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: muscle synergy, EMG-driven model calibration, synergy extrapolation, muscle excitation, stroke, muscle force

Citation: Ao D, Vega MM, Shourijeh MS, Patten C and Fregly BJ (2022) EMG-driven musculoskeletal model calibration with estimation of unmeasured muscle excitations via synergy extrapolation. Front. Bioeng. Biotechnol. 10:962959. doi: 10.3389/fbioe.2022.962959

Received: 07 June 2022; Accepted: 02 August 2022;
Published: 07 September 2022.

Edited by:

Li Li, Georgia Southern University, United States

Reviewed by:

Claudio Pizzolato, Griffith University, Australia
Zimi Sawacha, University of Padua, Italy

Copyright © 2022 Ao, Vega, Shourijeh, Patten and Fregly. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Benjamin J. Fregly, fregly@rice.edu

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.