Neural Data-Driven Musculoskeletal Modeling for Personalized Neurorehabilitation Technologies

Objectives: The development of neurorehabilitation technologies requires the profound understanding of the mechanisms underlying an individual's motor ability and impairment. A major factor limiting this understanding is the difficulty of bridging between events taking place at the neurophysiologic level (i.e., motor neuron firings) with those emerging at the musculoskeletal level (i.e. joint actuation), in vivo in the intact moving human. This review presents emerging model-based methodologies for filling this gap that are promising for developing clinically viable technologies. Methods: We provide a design overview of musculoskeletal modeling formulations driven by recordings of neuromuscular activity with a critical comparison to alternative model-free approaches in the context of neurorehabilitation technologies. We present advanced electromyography-based techniques for interfacing with the human nervous system and model-based techniques for translating the extracted neural information into estimates of motor function. Results: We introduce representative application areas where modeling is relevant for accessing neuromuscular variables that could not be measured experimentally. We then show how these variables are used for designing personalized rehabilitation interventions, biologically inspired limbs, and human-machine interfaces. Conclusion: The ability of using electrophysiological recordings to inform biomechanical models enables accessing a broader range of neuromechanical variables than analyzing electrophysiological data or movement data individually. This enables understanding the neuromechanical interplay underlying in vivo movement function, pathology, and robot-assisted motor recovery. Significance: Filling the gap between our understandings of movement neural and mechanical functions is central for addressing one of the major challenges in neurorehabilitation: personalizing current technologies and interventions to an individual's anatomy and impairment.


I. INTRODUCTION
H UMAN movement underlies interaction and coordination across the nervous, the muscular, and the skeletal systems [1]- [3] (see Fig. 1).Neural structures such as spinal motor neurons are ultimately recruited to produce mechanical forces [4], [5].Our understanding of human movement is still challenged by the difficulty of studying in vivo neural activity together with the resulting mechanical function elicited at the musculoskeletal Fig. 1.Synaptic inputs ultimately converge to pools of spinal motor neurons, whose spike trains trigger the innervated muscles.Spatiotemporally superposed muscle fiber action potentials result in EMGs with mechanical force produced after the electromechanical delay (EMD).Modeling formulations that reconstruct the subject-specific neuromusculoskeletal pathway can enhance our understanding of human movement and open up to novel HMIs.The short and stable muscle EMD compared to the larger and unstable corticomuscular delays make muscle electrophysiological recordings ideal for HMIs (see Section III-A).The EMD is sufficiently large (10 -80 ms) to enable models to predict the intended movement from EMGs and sufficiently consistent to be considered as a hard real-time deadline for the assistive device control.Meeting this deadline within an HMI would ensure the assistive device to be actuated synchronously with the user's muscle mechanical function.
level.Filling the gap between movement neural and mechanical functions is crucial for unraveling the subject-specific motor strategies that allow the composite neuromusculoskeletal system to operate across conditions such as motor tasks, training levels, impairments or pathologies.Understanding the mechanisms underlying an individual's neuromusculoskeletal function is central for addressing one of the major challenges in the field of neurorehabilitation technologies: personalizing the design of the assistive device, the human-machine interface (HMI), and the rehabilitation intervention to a specific patient's anatomy and impairment.
The neural processing underlying movement encompasses a number of stages at various levels of the nervous system, i.e., cortical and spinal.The resultant synaptic inputs ultimately converge to pools of spinal motor neurons, which perform the final neural processing step: transmitting the cumulative neural drive to the muscular level [3].The combination of a motor neuron and its innervated muscle fibers form a functional unit, called "the motor unit."This unit is the interface between the neural and the muscular levels and contributes to transduce neural inputs into mechanical outputs, i.e., motor neuron action potentials directly result in muscle fiber action potentials and, therefore, in mechanical force (see Fig. 1) [3].
The neural drive to muscles-the ensemble of action potentials of all active motor neurons-ultimately controls multiple musculotendon units (MTUs)-ensembles of muscle fibers and series elastic tendons.MTUs operate under neurophysiological constraints and transfer mechanical force to skeletal structures and degrees of freedom (DOFs) [6]- [9].The transformations involved are highly nonlinear (see Fig. 1).For the same neural command to a muscle, the force and moments transferred to the skeletal system can differ as a function of the muscle state, including fiber length, contraction velocity, or fatigue [2], [10].In this context, solving for the forces produced in the musculoskeletal system is an undetermined problem.There are more muscles than joint DOFs, so that an individual can use a variety of muscle recruitment strategies to produce the same joint moment and angle [2].
Establishing relationships between movement neural and mechanical functions is a complex long-standing problem.This challenge has been addressed using model-free machine learning approaches where all the transformations between given experimental variables, such as muscle electromyograms (EMG) and joint angles, are absorbed into a single transfer function or nonlinear mapping [11].An alternative approach is that of neuromusculoskeletal modeling, which characterizes each intermediate transformation by simulating the dynamic interplay between the nervous, the muscular, and the skeletal systems [6]- [8].
Neurorehabilitation interventions require clinically viable interfaces with the patient's nervous system and patient-specific mappings to the intended movement for effectively replacing or restoring the lost motor capacity.Muscle EMGs currently represent the only clinically viable solution for indirectly interfacing with the nervous system in patients, at least in the short to midterm.Invasive recording techniques from nerves and brain areas are currently limited mainly to academic scenarios [12]- [14].Due to the strong synaptic connection between a motor neuron and the innervated fibers, muscle fibers action potentials (mixed in the EMG signal) represent the biologically amplified version of motor neuron action potentials.Therefore, muscles can be seen as biological amplifiers of the synaptic input converging to the motor neuron pools (see Fig. 1).Surgical procedures such as targeted muscle reinnervation (TMR) [15] have further highlighted this concept.TMR is applied to individuals with proximal amputations (i.e., glenohumeral or transhumeral levels), where residual nerves that used to control distal muscles and DOFs are rerouted to residual proximal muscles, i.e., pectoral and dorsal muscles.In this scenario, the muscle is no longer used because of its mechanical function but solely for its ability of amplifying nerve signals associated to the function of the amputee's phantom limb.The amplified nerve signals are viably recorded via the reinnervated muscle EMGs [16], [17].
This paper aims to discuss clinically viable methods for extracting neural information from electrophysiological recordings and the development of subject-specific musculoskeletal modeling formulations that can be driven by such recordings.We will refer to this specific modeling class as to "neural data-driven musculoskeletal modeling."The review discusses how this modeling approach can be used for understanding the function of the human neuromusculoskeletal system in vivo in the intact moving human with emphasis on the applications that this can have in the context of neurorehabilitation and Fig. 2. Model-free approach first necessitates a learning step (A-1) where a mapping function is established between two given variables using for instance regression analysis.Only after the mapping is created (A-1), prediction of one variable as a function of the other (A-2) can be done.In the model-based approach, the analytical model can inherently predict "variable 2" as a function of "variable 1" (B-1).Prediction may be further improved by adjusting internal parameters as part of the parameter identification step (see B-2, Section IV).This is a closed-loop formulation where the model predictions (B-1) are continuously refined (or controlled) by model's internal parameters adjustments (B-2).The model-free approach does not rely on this closed-loop formulation, i.e., regression (A-1) and prediction (A-2) are two distinct and sequential steps.augmentation technologies.Despite the advances in the field, important challenges exist for translating EMG and model-based approaches to pathological populations and to the development of personalized and intuitive HMIs.Section II presents and compares model-free and model-based methodologies for linking between movement neural and mechanical function.Section III outlines how surface EMG can be recorded, processed, and used for interfacing with the human nervous system.Furthermore, it introduces techniques for recording and processing movement data and for establishing subject-specific musculoskeletal models.Section IV focuses on modeling techniques that use neural information extracted from EMG data to drive musculoskeletal plants.Section V describes how model-based methods can be employed to access internal neuromuscular variables that could not be viably measured experimentally (i.e., muscle force, joint moment, bone-to-bone compressive force, or joint stiffness) and use them for designing personalized assistive devices and HMIs.Sections VI and VII provide discussion and concluding remarks, respectively.

II. PREDICTING FUNCTION FROM NEURAL SIGNALS
Human movement is typically explored via experimental recordings reflecting neuromuscular and mechanical variables, e.g., electroencephalograms, EMGs, foot-ground reaction forces (GRFs), segmental body kinematics, or muscle ultrasound images [2], [3].However, the sole observation and analysis of experimental data does not allow establishing a direct cause-effect relationships across movement variables [18].Relations between observed variables can be generally identified using two approaches: the model-free and the model-based approach.Both approaches undergo learning and prediction steps, whose structure and differences are depicted in Fig. 2.

A. Model-Free Approach
Given experimental observations of the variables of interest, the model-free approach uses machine learning for approximat-ing a numerical function that maps between these variables [see Fig. 2(A)], i.e., from EMGs to joint angles.The mapping function can be approximated by combination of basis functions that can range from linear to polynomial and from exponential to Gaussian.Alternatively, artificial neural networks (ANNs) can be used for identifying mapping functions with no need for a specific basis of functions.
ANNs have been successfully used to approximate a variety of functional relationships underlying animal and human movement, i.e., from EMG to tendon forces in the cat soleus and gastrocnemious muscles [19], [20], to joint moments and joint angles, in the human lower extremity during locomotion [21].Linear and nonlinear regression has also been used to establish functional relationships between EMG signals and wrist joint angles [22]- [24] as well as hand kinematics [24], [25].In non-human primates, machine learning techniques have been frequently used to establish relationships between population of neurons experimentally recorded in the brain motor, premotor or parietal cortexes, and the final hand grasping action [26]- [28].
In the model-free approaches, functional relationships are established without explicit equations modeling the underlying neuromechanical processes.Therefore, the model-free approach is essentially a black box where all intermediate functional relationships between the observed experimental variables are not explicitly modeled but rather included into a macroscopic transfer function [see Fig. 2(A)].This also represents the main limitation of this approach.Functional relationships learned in a specific condition, or regime, may not generalize to novel conditions, as shown previously [11], [20], [29], [30].This may be related to the fact that a single macroscopic nonlinear transfer function may not be sufficient for capturing the complexity of the several intermediate nonlinear components existing between the observed variables.In this context, the more complex the regression map, the more it will overfit the data with resulting lack of generalization (i.e., classic overfitting problem).
Model-free methods do not enable direct understanding of the mechanisms underlying the observed variables.For instance, if EMGs are mapped into joint angles, a model-free approach does not reveal how limb movement emerges from individual's muscle coordination and function.Joint actuation from each spanning MTU can come from the series-elastic tendon or from the active muscle fibers [10].The mechanism of load sharing between muscles and how this results into joint kinematics is not revealed by the model-free approach.

B. Model-Based Approach
The model-based approach analytically defines each constitutive unit and functional relationships between the observed variables [see Figs. 1 and 2(B)].In the context of human movement, constitutive units can refer to neurophysiological and musculoskeletal components such as motor neurons, muscle fibers, series-elastic tendons, articular joints, or anatomical limbs.On the other hand, functional relationships define the way in which each constitutive unit converts input variables into output variables (see Fig. 1).Transformations, such as those from muscle force to joint acceleration, are described by analytical formulations with unknown parameters, i.e., the moment arm matrix (i.e., from muscle force to joint moment) or moment of inertia (i.e., from joint moment to joint acceleration).
Models of the composite neuromusculoskeletal system have been proposed in which the contribution of individual MTUs to joint actuation is resolved by using optimization.In these models, muscles are recruited according to the predefined optimization criteria such as minimizing the sum of squared muscle activations, so that the emerging muscle-actuated movement tracks experimental joint mechanics [31]- [33], or shows agreement with joint mechanics normative values [34], with additional control objectives such as minimal metabolic cost of transport [35].Alternative methods have been proposed in which the contribution of individual MTUs to joint actuation and whole-body locomotion is resolved by solely modeling muscle and spinal reflexes [9], [36] and by optimizing muscle-tendon control parameters to minimize the metabolic cost at different locomotion conditions [37].Models based on networks of spiking motor, sensory, and interneurons have been recently used for the control of simple single-joint musculoskeletal models during postural tasks [38], [39].Models of the animal and insect nervous system have been proposed that simulate the role of spinal central pattern generating networks in the control of muscle recruitment rhythm and in the coordination of left-right movements [40]- [42].Comprehensive modeling formulations have been proposed that simulate the interaction between central pattern generating networks in the nervous systems with local feedback from sensory neurons encoding kinematic and dynamic information generated at the musculoskeletal level [43]- [45].
These are fully predictive modeling formulations that simulate how neural structures interact with the musculoskeletal plant and enable predicting emerging motor behaviors in unknown/unmeasured scenarios.Although current neuromusculoskeletal models provide a valuable position for the computational investigation of neuromuscular strategies [46], [47], an individual's neuromusculoskeletal function (i.e., neuronal interaction, motor unit recruitment, musculoskeletal force production) is highly variable across motor tasks [48]- [50], pathology [51]- [53], or training [54], [55].Differences in neuromusculoskeletal function occur across individuals [56] but even within the same person [49], [50].This limits current models' ability of producing valid and realistic neuromechanical predictions that can account for conditions including neuromuscular deficits, orthopedic impairment or an individual's neuromusculoskeletal response to neurorehabilitation interventions or external assistive devices [57], [58].
Muscle electrophysiological recordings (i.e., EMG) in conjunction with advanced processing techniques can be used to gain experimental access to the neural information underlying an individual's movement (see Section III-A) [1].This enables deriving experimental estimates of subject-specific and contextdependent neuromuscular strategies with no direct need for creating numerical models of these.In this scenario, EMG estimates of neuromuscular excitations have been used to drive neuromusculoskeletal models during a variety of dynamic motor tasks as well as neuromuscular [59], [60] and orthopedic conditions [51], [61]- [63] and predict a range of neuromuscular variables including: muscle forces [51], [64], [65], joint moments [6], [51], Fig. 3. Movement data are collected from the moving subject (see A, Section III).Retroreflective markers are placed on the human segments and their 3-D positions are tracked using optical motion capture cameras with sampling rates typically between 50 and 250 Hz depending on the dynamic of the motor task [2].GRFs are recorded (2 KHz) using in-ground force plates.Marker trajectories and GRFs data are low-pass filtered with the same zero-phase filter and cut off frequency in order to preserve dynamic consistency across data.Raw EMG signals are processed (B) to extract neural excitation estimates including: linear envelopes (blue curves in B), lower-dimensional synergy modules (green Gaussian primitives and weighting coefficient in B), or the EMG constituent motor unit discharge events (red impulses in B).Reflective marker trajectories during static poses are used to scale a musculoskeletal template (C).The scaled model is used to extract joint angles (D) and moments (E) using IK and ID or direct sensor recordings.Joint angles and neural excitation estimates are used to drive the neural data-driven model (F).The model's parameters are calibrated to an individual (G) using identification methods (see Section IV).Validation is done (G) on the model's ability of blindly predicting experimental data (i.e., joint moments, compressive forces, angles, etc.) during novel trials.
In general, the model-based approach is directly limited by uncertainties arising from the fact that we cannot always measure internal physiological parameters or rely on the existence of "proven laws" that accurately describe the dynamics of the constitutive units and functional transformations involved.Often, such units and transformations are built upon hypotheses, rather than laws.In this context, model estimation-detection theory is used to identify the unknown model parameters [see Section IV, Fig. 2(B)] [11].

III. DATA COLLECTION
Subject-specific neural data-driven musculoskeletal models can be established from experimental data that typically include EMG, joint kinematics, and joint dynamics [see Fig. 3(A)] [76].

A. Extracting Neural Information From the Surface EMG
Surface EMGs contain information on the ensemble of spike trains discharged by motor units (see Fig. 1).Direct or indirect estimates of this information will be termed in the reminder of this paper as "neural excitation" [see Fig. 3(B)].This section covers methods for extracting neural excitations from multimuscle EMGs, from low-dimensional approximations of these, and from the EMG constituent motor unit action potentials.
Estimates of neural excitations are conventionally derived from the intensity (i.e., amplitude) of surface EMG [77].Raw EMG signals need to be demodulated to recover amplitude information proportionally to the neural drive to the muscle [78].For this purpose, raw EMGs are sampled (typically at 1000 Hz) and high-pass filtered to remove low-frequency noise.The cutting frequencies can largely vary between 10 and 250 Hz with higher cut off frequencies being suggested to better represent the contribution of distant motor units to the EMG amplitude [79].The resulting signal is full-wave rectified and low-pass filtered (2 -6 Hz) [7], [79], [80].Low-pass filtering is performed because muscle force operates with lower dynamics (i.e., it has a smaller bandwidth) than the high-pass filtered and rectified EMG signal.This is a common feature of mechanical motors and in biological actuators, such as muscles, it reflects the dynamics of electrochemical transformations including calcium dynamic, limits in the propagation velocity of fiber action potentials, and muscle-tendon viscoelasticity [2], [3], [8].Therefore, low-pass filtering enables better correlating the resulting EMG-linear envelope to muscle force.This aspect was further confirmed at the level of motor units showing that the low-frequency oscillations in the motor output from the spinal cord best correlates with the overall force produced by the muscle [81].EMG-linear envelopes are then normalized with respect to the peak-processed values obtained from maximal voluntary contractions or across entire sets of recorded dynamic trials [see Fig. 3(B)].Obtaining true EMG maxima from maximal voluntary efforts is a critical step and a source of uncertainties [3].
An alternative to classic filtering for EMG amplitude estimation is the Bayesian filtering, such as Kalman filtering [82], particle filtering [83] or Gaussian process regression models [84].Recursive Bayesian filters can be used to estimate in real time the state of a dynamic system, i.e., surface EMG amplitude (state) in a specific muscle (system) [85].These filters are based on probabilistic prediction and observation models and are learned from a set of training data.In the learning step, the probability distribution used to solve for the Bayes theorem is identified and used to subsequently estimate the EMG amplitude.The advantage of these filters is that the cut off frequency does not need to be explicitly defined a priori, but is dynamically determined based on the models learned on the training data.This enables a tradeoff between smoothness of estimates and the ability of predicting rapid amplitude changes in the extracted envelope [84], [85].EMG amplitude estimation via classic or Bayesian filtering is the most common approach for extracting neural excitation estimates and for generating control commands for subsequent modeling stages (see Section IV, Fig. 3).These methods assume that the estimated EMG amplitude is linearly proportional to the actual motor unit firings.This assumption is, however, not fully met because of several sources of uncertainty including electrode placement, detection volume, and motor unit action potential superimposition.
Recently, other approaches have emerged that enable extracting EMG features that better reflect the muscle neural excitation, either at the macroscale (i.e., muscle synergies) or at the microscale (i.e., motor units).In this scenario, unsupervised learning techniques are used to extract linear [86]- [89] or nonlinear [90], [91] approximations of the structure and statistical distribution underlying the patterns of multimuscle EMGs [see Fig. 3(B)].Linear approximations of low dimensionality give the possibility of capturing the elementary modules (or synergies) of multimuscle coexcitation that are descriptive of a specific regime or motor condition.Increasing experimental evidence is supporting the hypotheses that EMG-extracted synergies well represent the neural coordinative structures, explicitly encoded in the human nervous system, that are used to reduce the computational burden in the neuromuscular control of complex motor tasks [92]- [94].Synergies can be extracted using methods such as principal component analysis [95], independent component analysis [96], or non-negative factorization [86], [97].The extracted synergies can well describe the statistical distribution of muscle excitations and account for their variability.Importantly, they can be used to synthetize the process of muscle recruitment under specific regimes and enable driving musculoskeletal models with a small number, or without, experimental EMG data [see Fig. 5(B)] [98]- [100].
Nonlinear transformations of EMG signals give the possibility of identifying the underlying processes responsible for the generation of the surface EMG.These methods include blindsource separation [90], [101] and Bayesian decomposition [91] of high-density multichannel EMG signals.Blind-source separation is now an established powerful approach for decomposing the EMG signal into the constituent motor unit activity and for resolving the inherent problem of motor unit action potential superposition [see Fig. 3(B)] [101].This relies on a mathematical convolution model with additive noise that describes the mixing process of motor unit action potentials and spike trains in the generation of the interference EMG signal [102].Once the mixing model is established, the inverse model can be applied to experimentally recorded EMGs, thus enabling mapping from experimental (mixed) EMG data to (unmixed) motor unit spike trains [see Fig. 5(B)] and action potentials.In this scenario, the possibility of recording high-density EMG signals, typically using two-dimensional (2-D) grids of densely located electrodes, substantially increases the number of detectable motor units with respect to recordings based on nondense electrodes [101], [103], [104].
High-density EMG and decomposition provide experimental access to motor unit spiking events, which corresponds to those produced by spinal motor neurons [1].This is central for investigating motor control strategies and motor unit morphological and functional properties, as well as their adaptation to conditions such as pathology, fatigue, pain or exercise.Importantly, this directly addresses uncertainties in current EMG amplitude estimation techniques with profound implications for realizing myoelectric HMIs.
The signal processing methods presented in this section provide experimental access to neural information in vivo in the intact moving human.This approach is particularly powerful when combined with musculoskeletal modeling formulations (i.e., neural data-driven modeling, Section IV), thus enabling personalized neurorehabilitation treatments and HMIs (see Section V).

B. Kinematics and Kinetics of the Human Skeletal System
Stereophotogrammetry in conjunction with models of the human skeletal mechanics [HSM, Fig. 3(D), (E)] represent the gold-standard solution for deriving the mechanics of anatomical joints and limbs [2].Retroreflective markers and GRFs are recorded and processed as described in Fig. 3.The HSM model numerically represents anatomical segments (i.e., limbs) interconnected by joints, influenced by forces (i.e., GRFs), and subject to constraints (i.e., joint range of motion) [105].Computationally efficient implementations treat anatomical segments as rigid bodies that do not undergo deformation when force is exerted on them [106]- [108].Subject-specific HSM models are created from medical imaging data of bone and muscle surfaces, such as magnetic resonance imaging (MRI) or computed tomography [109].Freely available software packages such as the Musculoskeletal Atlas Project Client1 [110] and NMS Builder2 [111] enable rapidly generating detailed 3-D models using image fusion, segmentation, and mesh generation methods.
However, subject-specific HSM models are most commonly created using HSM templates further scaled to match the individual's anthropometry [see Fig. 3(C).These templates are obtained from data averaged across a large number of cadaveric specimens.They have the disadvantage of not accounting for subject-specific bone and muscle geometry, which might be particularly important to capture features of pathological populations.However, they have the advantage of enabling the application of musculoskeletal modeling to any subject whether or not medical imaging data are available to begin with.This is important in the context of neurorehabilitation technologies.Scaling the HSM template to an individual's anthropometry is done based on the marker trajectories recorded during a subject's static standing pose [see Fig. 3(C)] [112].In this process, the mass properties of the subject are also proportionally adjusted in order to match the subject's total measured mass.The scaled HSM model is then used to reconstruct the 3-D joint angles and moments based on the marker trajectories and GRFs recorded during dynamic motor tasks by solving inverse kinematics [IK, Fig. 3(D)] and inverse dynamics [ID, Fig. 3(E)], respectively.The steps described in this section can be performed using the modeling software package OpenSim3 [112].Additional software packages include AnyBody 4 and Biomechanics of Bodies. 5he methods so far presented rely on nonportable in-ground force plates and multicamera systems.These are not easily applicable in the context of neurorehabilitation or augmentation technologies, which often require fully wearable solutions.In this context, the mechanics of human skeletal joints can be measured using wearable angle and force sensors located in the correspondence of a specific joint center of rotation [see Fig. 3(A)].However, human joints do not have, in general, fixed centers of rotation nor orthogonal axes [2], [113], [114].Furthermore, the morphology of human segments largely varies across individuals [2], [115].These factors limit the ability of current technology to directly measure accurate joint dynamics [116].The combination of wearable sensors such as inertial measurement units and insole pressure sensors with model-based IK and ID techniques may provide fully portable solutions capable of extracting accurate joint dynamics estimates with no substantial loss of accuracy with respect to gold-standard methods [117]- [119].

C. Kinematics of the Human Muscular System
From the imaging data used to create skeletal bone and joint representations (see Section III-C), the origin and insertion points of selected MTUs can also be segmented [109].The shortest line segment connecting the MTU origin and insertion points via muscle-to-bone wrapping points is typically used to represent the geometry of fibers in series with tendons as a unidimensional wire-like path [109].This has been shown to accurately represent the kinematics (i.e., length, velocity, and moment arms) of MTUs, both in the upper [120] and lower [121] extremities.Alternatively, finite element modeling can be used to create 3-D volumetric muscle representation from muscle imaging data [122]- [124].However, these methods are typically associated with large computational costs that increase with the volumetric surface to be represented.
The MTU kinematics derived from 3-D musculoskeletal models can then be synthetized into simpler surrogate models using computationally inexpensive structures such as multidimensional cubic B-splines [125], or polynomials [126] that operate as a function of joint angles.Spline-based estimates of MTU kinematics were recently used to simulate musculoskeletal kinematics in real time on embedded systems with limited computing capacities with no need for simulating complex 3-D skeletal and muscular geometries [64], [127].

IV. NEURAL DATA-DRIVEN MUSCULOSKELETAL MODELING
Current sensing and modeling technologies provide solutions for recording muscle EMGs (see Section III-A), joint kinematics and kinetics (see Section III-B), and MTU kinematics (see Section III-C).These sensing and modeling steps can now be performed in real time using fully wearable technology, which is an important prerequisite for neurorehabilitation technologies.However, in order to determine the force produced by the musculoskeletal system, it is necessary to account for how muscles are neurally recruited and how they transfer force around multiple joints.To address this problem, EMG estimates of neural information [see Section III-A, Fig. 3(B)] and joint kinematics can be used within a computational framework that includes models of muscle neural activation, kinematics, contraction dynamics, and joint mechanics (see Fig. 4).Open-source toolboxes such as CEINMS 6 [128] or EMG-FE7 [129] are currently avail-  [125].The MTU dynamics component (C) solves for the dynamic equilibrium between muscle fibers and series elastic tendons in the production of MTU force [7], [176].It employs a Hill-type muscle model informed by MTU length and neural activations from the previous two components.The joint dynamics component (D) transfers MTU forces to the skeletal joint level using MTU moment arms.In the offline calibration component (E) initial nominal parameters are repeatedly refined, as part of a least-squares optimization procedure, so that the mismatch between model's predicted joint function and the experimentally recorded joint function is minimized [7], [8], [65].
able that may boost the adoption and development of these techniques.

B. Open-Loop and Closed-Loop Formulations
The calibrated model can be operated as a part of an openloop formulation (see Fig. 4) [7], [65], [98].This is a "predictive formulation," where neural excitation estimates directly drive MTUs to compute the resulting joint mechanical function blindly, i.e., with no tracking mechanism that accounts for prediction errors [see Fig. 5(C)].This calibrated open-loop formulation can be generalized into a closed-loop formulation, also regarded as a "tracking formulation" [6].This was recently proposed to account for uncertainties in neural excitation estimated as EMG-linear envelopes including cross-talk, processing artifacts, and the inability of recording surface EMG from deep muscles.These are factors that may limit the open-loop model ability of matching experimental multijoint moments precisely [see Fig. 5(C)] [6].To account for them, the calibrated openloop modeling formulation (see Fig. 4) is coupled with a static optimization block [31].The resulting model allows extracting neural excitations by balancing information from experimental EMG-linear envelopes and those derived from static optimization solutions.
It was recently shown that a minimal adjustment of experimental EMG-linear envelopes (i.e., root-mean-squared error < 5%), in conjunction with the static optimization solution for the hip flexing iliopsoas MTU, was sufficient to precisely track experimental multi-DOF joint moments during walking and running, with substantial improvement with respect to open-loop formulations.Importantly, this also allowed characterizing the dynamics of joints such as the hip for which the major flexors (i.e., iliacus and psoas MTU) are deeply located and inaccessible via surface EMG [see Fig. 5(A)] [6].This provided for the first time simulations that were consistent with muscle electrophysiological activity [see Fig. 5(A)] and with movement mechanics [see Fig. 5(C)].

V. APPLICATIONS TO NEUROREHABILITATION TECHNOLOGIES
In the reminder of this section, we will present examples of application areas where neuromusculoskeletal modeling is highly relevant for establishing HMIs and for characterizing the dynamics of human-machine interaction.This section will present applications of neural data-driven modeling formulations where neural excitations are estimated as EMG-linear envelopes.This special case will be referred to as "EMG-driven modeling" [7], [8], [65].

A. Orthoses
Orthoses are used to restore or augment musculoskeletal function in impaired or healthy individuals.The main limitation in commercial ambulatory powered systems (e.g., ReWalk Robotics, Ltd., Israel) is that the available HMIs permit little or no user's active participation in the device control, thus limiting neuromuscular plasticity and the subsequent motor function recovery.Furthermore, there is currently no commercially available orthosis that can preserve or decrease the metabolic cost of locomotion in healthy subjects.This suggests that current solutions may trigger abnormal and suboptimal musculoskeletal function [137].
A number of model-based solutions have been proposed to address these limitations.The works in [138] and [139] presented powered orthoses for upper and lower limbs that modulate the support to the user based on the user's predicted strength, i.e., the ability of producing joint moment.This is estimated in real time directly from EMG-linear envelopes and wearable joint sensors using an open-loop EMG-driven modeling formulation (see Section IV-B, Fig. 4).This allows establishing controllers that actuate the orthosis proportionally to the EMG-predicted user's joint moment capacity.The predicted moments can be augmented by a support ratio to define the desired level of joint support that the orthosis will provide to the user (see Fig. 6) Fig. 6.Model-based HMI scheme for lower limb powered orthoses and variable-stiffness prostheses.(Orthoses) The user's joint moment generating capacity is estimated from the affected leg EMGs, with the powered orthosis being actuated proportionally and with support levels varying as a function of the "support ratio" coefficient [138].(Prostheses) Physiological joint stiffness profiles are estimated off-line from healthy subject populations or on-line from unaffected leg EMG recordings.The prosthesis is then controlled to reproduce stiffness templates [155] (see Section V). [138].A simplified modeling formulation was used in [140] to predict in real time the user's knee joint stiffness as a function of EMGs and modulate the impedance of a compliant knee joint orthosis accordingly.These examples based on EMG-driven modeling allow establishing a symbiotic HMI, where the powered orthoses dynamically adapt to the user's motor ability (see Fig. 6).In this scenario, generating the device control command within the muscle electrophysiological delay would enable the worn device to naturally respond to the user's neuromuscular function (see Fig. 1).This delay can range from shorter intervals (i.e., ∼10 -30 ms) typically in lower limb muscles [141] to larger ones (i.e., ∼ 80 -100 ms) typically in upper limb muscles [142].
An unpowered ankle joint orthosis with a spring mechanism acting in parallel to the calf muscles was presented in [143] that reduced for the first time the metabolic cost of human walking by ∼ 7%, something remarkable given the optimality of human locomotion.Interestingly, metabolic energy savings could be achieved only for specific intermediate orthosis spring stiffness, i.e., too stiff or too compliant springs did not provide substantial benefit.In [144], musculoskeletal modeling was used to directly link musculoskeletal mechanics and energy consumption when hopping with a similar orthosis.Results showed that suboptimal spring stiffness levels in the elastic ankle orthosis caused unfavorable soleus muscle-tendon mechanics, i.e., increased soleus contraction velocity around unfavorable ranges of fiber length.The suboptimal soleus mechanics determined the increased metabolic cost of motion [143], [145].Further studies [37], [146] used musculoskeletal modeling to show that metabolic efficiency in locomotion is tightly related to plantarflexor muscles operating with minimal contraction velocity.Although not used for establishing an HMI, modeling was used in these studies to understand musculoskeletal function underlying efficient locomotion and efficient human-machine interaction.
The ability of modeling and simulating the interaction of humans with devices connected in parallel with their limb is central for deploying wearable technologies that can effectively enhance an individual's mobility.Fig. 7(A) shows an example where musculoskeletal modeling was employed to study human-machine interaction in a subject with severe quadriceps Fig. 7. (Left) (A) Musculoskeletal model of a subject with severe quadriceps tendon tear and with parallel passive knee orthosis.The subject's quadriceps force was set to zero in the model.The orthosis was modeled using knee joint force constraints acting in parallel to residual muscle force.This characterized the orthosis joint stiffness response to joint bending and enabled keeping the subject's knee joint extended during the gait stance phase.(Right) (A) Measured and predicted joint load carried by the orthosis during locomotion (see Section V-A).(Left) (B) Musculoskeletal model of a transfemoral amputee with series prosthesis.Tibiofemoral compressive forces from the intact contralateral knee were estimated using OpenSim joint reaction analysis tool [112].(Right) (B) Compressive force expressed as body weight multiples (xBW) predicted for the amputee (see Section V-B) and from one healthy control individual as well as in vivo measurements from 10 individuals with force-measuring knee implants (reference).In both A and B, OpenSim was used to build a musculoskeletal model that accounted for the subject's anthropometry as well as that of the worn device [112].A closed-loop EMG-informed modeling formulation (see Section III-B) [6], [46] was used to realize a muscle-driven simulation that tracked experimental data including experimental EMG, joint moments and joint angles.Simulations are provided here as representative applicative examples rather than complete validation results.Movement data were recorded and kindly made available to the authors by the OttoBock Health Care Biomechanics Laboratory (Göttingen, Germany).tendon tear.The subject walked wearing a passive knee-anklefoot orthosis (OttoBock, Germany) that locks the knee joint during the stance phase, preventing the subject's knee to collapse.The orthosis-assisted locomotion was simulated using Open-Sim [112] and a closed-loop EMG-informed modeling formulation (see Section IV-B) [6].Fig. 7(A) shows how this modeling framework enables us to accurately predict the interaction forces between the orthosis and the subject's leg, i.e., the load carried by the passive orthosis knee joint.Experimental interaction forces were measured by a wearable force-measuring sensor (Otto-Bock, Germany) mounted on the orthosis knee joint.In [147], an optimization framework was developed using OpenSim [112] to explore the optimal design of simulated wearable robotic devices that maximally enhanced long jump performance.

B. Prostheses
Prostheses are used to replace anatomically or functionally missing biological limbs.The main limit in current prosthetic systems is the control scheme.Current control schemes for commercial myoelectric upper limb prosthesis are still unintuitive and provide limited regain of lost functionality [148].This is the main reason for the high abandonment rate of current commercial solutions [149].While model-based approaches to the myoelectric control of upper limb prostheses are yet to be proposed, a number of model-free solutions (see Section II-A) has been presented in the past decades based on classification and regression (see Section II-A) [148].Although these methods enable more natural (i.e., proportional and simultaneous) control of multiple joint DOFs than current commercial solutions, they support limited arm and hand functions and are validated in specific laboratory settings [14].
Lower limb prosthesis technology has made major advances in the past decades substantially improving amputees' mobility and fall risk prevention.Latest commercial solutions (i.e., Ottobock Genium and Össur Rheo Knee) employ microprocessorcontroller joint dampers that improve gait performance over purely passive prostheses [150].However, current prostheses are still inferior to their biological counterparts in that they neither generate positive net work over the gait cycle, nor modulate joint torques and viscoelastic forces as naturally as humans do [151].As a consequence, lower limb amputees who reply on prostheses suffer from increased energy consumption during gait and are challenged when locomoting across different terrains or transitioning across motor tasks or locomotion speeds.Furthermore, the functional discrepancy between the artificial and the biological leg results in asymmetric gait with abnormal loads on the healthy leg, which may trigger degenerative bone disorders such as osteoarthritis (OA) [152]- [154].Model-based solutions for lower limb prostheses have been recently proposed that are addressing current open challenges.
In [71], a simplified open-loop EMG-driven formulation (see Section IV-B) was used to extract values of muscle short-range stiffness and how it projects onto the knee joint across locomotion tasks performed by healthy intact individuals.These were subsequently used in [155] to generate knee joint stiffness template values for designing a biomimetic variable stiffness transfemoral prosthesis (see Fig. 6).The resulting prosthesis displayed task-adaptive behavior, i.e., it was capable of modulating joint stiffness (derived from muscles short-range stiffness), moment, and velocity across different motor tasks.In [156] and [157], a transtibial ankle prosthesis controller was developed based on a reflex-based neuromuscular model of the ankle joint muscles.The model simulates local reflex loops (i.e., muscle stretch reflex) on the basis of the prosthesis joint kinematics and interaction with the environment.The simulated reflexes recruit the virtual muscles and the resulting simulated muscle forces are used to control the prosthesis joint in a human-like manner.Results showed that the prosthesis displayed speedadaptive behavior, i.e., it could modulate net joint work across a range of locomotion speeds.A similar reflex-based neuromuscular modeling approach was used in [158] to control a transfemoral knee-ankle-powered prosthesis with improved functionality for recovering balance from unexpected external disturbances.Preliminary results in simulation and from healthy individuals showed the model-based prosthesis control allowed for walking on rougher terrains and for recovering from larger trip disturbances than conventional impedance controlled prostheses.
Likewise for the orthosis scenario (see Section V-A), understanding how the amputee's musculoskeletal system responds to an artificial leg connected in series with their residual limbs is crucial to deploy prostheses that can come closer to biological limb performances.Fig. 7(B) shows an example where musculoskeletal modeling was used to understand humanmachine interaction in the context of an amputee walking using a commercially available transfemoral prosthesis (Ottobock, Germany).This was simulated using OpenSim [112] and a closed-loop EMG-informed modeling formulation (see Section IV-B) [6].Fig. 7(B) directly compares the tibiofemoral forces extracted from the amputee with those measured in vivo from ten subjects that were implanted with a force-measuring knee replacement [133].The in vivo force data were taken from the OrthoLoad database. 8Simulations suggest how amputees' contralateral knees may carry larger loads than healthy individuals' knees, throughout each gait cycle.In this example, musculoskeletal modeling provides the ability of predicting internal variables, which are difficult to measure experimentally (i.e., tibiofemoral forces), as well as an estimate of how the human musculoskeletal system may interact with and respond to artificial limbs.Although further validation is required, this approach offers invaluable opportunities for designing future prostheses that can effectively restore contralateral limb function.

C. Gait Retraining
Gait retraining is a noninvasive treatment strategy for altering and controlling target gait variables associated with the progression of orthopedic or neurological conditions [159], [160].It has recently received attention, especially in the treatment of orthopedic conditions, such as OA.In this application, it is proving to be an effective treatment for correcting gait alterations that often persist after surgical interventions [161] or as an alternative to such interventions, i.e., total joint replacement [160].Current gait retraining methods for knee OA rely on the use of simple biomechanical models for calculating the external knee adduction moment (KAM) as a target variable to control during the gait retraining interventions [160].Decreasing the early stance peak KAM has been reported to also decrease pain, disease progression, and disease severity in OA patients [154], [162], [163].Recent studies have explored real-time visual and vibrotactile feedback to enable subjects to relearn their gait with reductions in KAM that ranged from 7% to 48% [160], [162], [164].Similar procedures are currently being developed and tested for treating hip OA [165].
Gait retraining has been also prescribed to patients with forcemeasuring knee implants, where the target variable to alter (i.e., to minimize) was the in vivo tibiofemoral contact force [166].Retrained gaits with minimal in vivo tibiofemoral contact forces may be more effective than gaits with minimal KAM peaks to the treatment of OA condition [68] because OA progression is directly related to tibiofemoral contact forces and only indirectly to KAM peaks [68].The availability of EMG-driven models that can predict accurate estimates of tibiofemoral forces [68] in real time [64], [167] will offer in the immediate future the possibility of performing joint contact force-based gait retraining to any subject, even in those who are not fitted with force-measuring implants.This is also expected to enable broader injury prevention applications where estimates of internal neuromuscular variables (i.e., muscle-tendon strain) are continuously monitored and controlled to stay within safe limits.

VI. DISCUSSION
We have presented an overview of clinically viable methods for interfacing with the human nervous system in vivo and for modeling the resulting function of the musculoskeletal system.In the context of neurorehabilitation technologies, the challenge is to build neuromusculoskeletal models that are both computationally efficient (i.e., for real-time HMIs) and physiologically correct (i.e., for characterizing the subject).
We first introduced the concept of muscles as biological amplifiers of the spinal cord neural output (see Section I).In this context, we presented methods for extracting movement neural information from muscle EMGs (see Section III-A, Fig. 3).Although implanted nerve or cortical electrodes provide a rich source of neural information, these recording modalities are mainly limited to research scenarios [14], [168] and can be excluded for interventions of limited duration, such as robot-aided rehabilitation.The surface (or intramuscular [13], [169]) EMG currently represents the only viable solution that can be directly applied to existing clinical settings [14].This motivated this review's focus on techniques for extracting neural information from surface EMG data.
We then showed how neural features extracted from muscle EMGs can be further processed using neural data-driven modelfree or model-based approaches for the purpose of decoding an individual's motor intention (see Section II and Fig. 2).In the context of neurorehabilitation, these paradigms offer the unique opportunity of interfacing with the patient's nervous system in vivo and predicting the intended movement for best replacing/restoring the impaired motor ability.The model-free approach has the advantage of not requiring explicit mathematical models that mimic the nature of the intermediate transformations, i.e., from muscle EMG to limb movement.However, it does not allow understanding the mechanisms underlying the learned relationship as well as how they vary with pathology or changing muscle properties.Moreover, functional relationships learned in a specific condition, may not easily generalize to novel conditions [11].This limits the ability of understanding movement and, therefore, delivering personalized neurorehabilitation treatments and technologies.
Neural data-driven model-based approaches overcome this limitation through explicit analytical formulations that link electrical (neuromuscular) signals into muscular and articular mechanical forces (see Figs. 1 and 5).The main advantage of combining electrophysiological muscle recordings with comprehensive biomechanical models is that it enables accessing a larger spectrum of neuromechanical variables than analyzing electrophysiological data or movement data individually.That is because it enables establishing functional relationships in human movement such as assigning a direct "mechanical meaning" to the recorded electrophysiological data (see Section II and Figs. 1 and 5).Furthermore, it allows predicting the musculoskeletal response to any recordable neuromuscular behaviors with no direct need for creating numerical models of these.This is central for understanding and diagnosing neuromuscular deficits as well as for developing personalized neurorehabilitation treatments and intuitive HMIs.
Alternatively, fully predictive modeling formulations that do not rely on electrophysiological recordings offer the opportunity to perform broader "what-if" analyses by changing the simulated neural commands and by observing the emerging motor behavior [58].Challenges to this approach include the development and validation of fully predictive models of the human nervous system that can adapt over short-time scales (i.e., by triggering neuromuscular reflexes to external stimuli) or over long-time scales (i.e., for learning new motor programs by reorganizing spinal and brain structures) so that they can be used to predict valid responses to truly unknown and unmeasured scenarios (i.e., neuromuscular deficits) [58], [170].
Both neural data-driven and fully predictive modeling formulations have been exploited in the context of neurorehabilitation technologies (see Section V) for establishing myoelectric HMIs (see Fig. 6), for synthetizing biological limb function (i.e., joint stiffness modulation) into artificial limbs, or for characterizing human-machine interaction (i.e., Fig. 7).Section V outlined representative examples of myoelectric HMIs that employed EMG-driven modeling for controlling powered orthoses and prostheses (see Section V, Fig. 6).For this purpose, we have previously presented EMG-driven modeling formulations that can operate within physiological electromechanical delays (EMDs) [141] by employing computationally efficient models of the musculotendon kinematics [125] and dynamics [64] directly on embedded systems.Especially for powered orthoses, HMIs that actuate the wearable device solely on the basis of the detection of externally measurable forces (i.e., external joint moments or limb-orthosis interaction force) could not provide support until the user has produced detectable interaction force or movement [171].On the other hand, the ability of predicting internal musculoskeletal forces from EMG data enables predicting the user's movement intention before the movement actually takes place, thus enabling support even in individuals with reduced motor abilities but with detectable electrophysiological activity.
Another important aspect to consider for myoelectric modelbased HMIs is that of reducing the number of needed EMG electrodes.The accurate model-based prediction of forces in the lower extremity joints would require ∼ 15 EMG channels per leg [65].A wearable device driven by a large number of electrodes would become susceptible to sensor noise and artifacts, especially if used during dynamic limb motion.This may be the case for the myoelectric control of ambulatory powered orthoses and may be less of a problem in the control of upper limb prostheses where residual muscles operate under quasi-isometric conditions due to the amputation.A number of solutions exist for reducing sensitivity to sensory noise.Muscle coordination can be synthetized by employing the theory of muscle synergies [87].Muscle synergies-based models have been recently proposed for predicting muscle excitation [99] and the resulting musculoskeletal forces [98] [see Fig. 5(C)] during large repertoires of locomotion tasks with no loss of accuracy with respect to using experimental EMGs.This enables predicting muscle coordination as a function of joint kinematics and estimates of the gait cycle percentage, by assuming that muscles are recruited according to the predefined synergistic schemes [98].Alternatively, models of local neuromuscular reflex loops can be used to generate synthetic activation signals for the major lower limb muscles [172].Less computationally efficient solutions can rely on closed-loop EMG-informed simulations [see Section IV-B, Fig. 5(A)], where static optimization is solved on a frame-byframe basis for deriving force estimates for muscles with no EMG data available [6].This can be done in real time if nonembedded processing units are used.These approaches could be used to completely replace the need for EMG recordings or for decreasing the number of EMG sensors.
This review also showed the use of neuromusculoskeletal modeling for understanding cause-effect relationships in human-machine interaction (see Section V, Fig. 7).In this scenario, modeling and simulation are especially valuable for finding optimal device design that is personalized to an individual, i.e., finding the prosthetic limb properties that result in minimal contralateral knee joint loads [see Fig. 7(B)] or the orthosis minimal sizing properties for carrying the required load during stance [see Fig. 7(A)].This allows addressing the limitations of the traditional trial-and-error approach, where orthopedic devices are designed irrespectively of the user's musculoskeletal properties and optimized based on limited experimental data.
As discussed in Section II-B, the model-based approach is limited by the impossibility of measuring (and thus modeling) internal physiological parameters and variables, such as muscle optimal fiber, physiological cross-sectional areas or tendon slack length.This limit is addressed by using model-estimation theory as outlined in Section IV (see Fig. 4).In the context of patients, the uncertainties to account for may increase due to factors related to the orthopedic or neurological condition, which pose additional modeling challenges.In these scenarios, highfidelity patient-specific models describing pathological muscles and bone structures may be first created based on MRI imaging data and finite element modeling (see Section III.C-D) [173].These models could then be synthetized into simpler surrogate models that describe the complex musculoskeletal dynamics using computationally efficient structures such as multidimensional splines [125], polynomials [126] or ANNs [174] (see Section III.D).In this scenario, the ability of merging the modelfree and model-based approaches together will prove valuable, especially to compensate for patient's features that cannot be explicitly modeled.
It is important to critically assess the modeling complexity that is truly needed to address a specific problem.The examples of this review (see Section V, Fig. 7) showed how relatively simple models could be used to characterize patients with amputations or muscle weakness wearing assistive devices such as orthoses and prostheses.Furthermore, EMG-driven modeling formulations based on simple musculoskeletal geometries were successfully used in the past to characterize the neuromuscular mechanisms underlying individuals affected by orthopedic of neuromuscular conditions including anterior cruciate ligament rupture [61], patellofemoral pain [51], OA [62], [63], and stroke [59], [60].
Moreover, the availability of "big data" of human movement will help to further address the modeling limits.In a world where wearable technology (i.e., textile sensors, smart phones and watches, and orthotic devices) and cloud data are becoming increasingly pervasive, the availability of human movement data will reach larger scales than ever [175].Frameworks that can integrate modeling paradigms to process a variety of data types (i.e., high-fidelity data from laboratories or low-fidelity data from wearable embedded sensors) will prove critical to characterize movement function and pathology with statistical inference power and enhanced predictive capacity.This will enable making increasingly accurate predictions in scenarios where data are unavailable, i.e., predicting the progression of neuromuscular or orthopedic diseases, the response to a surgery or to physical training.This approach is especially promising for understanding and forecasting history-dependent and contextdependent neuromusculoskeletal function (i.e., fatigue, training, and impairment).

VII. CONCLUSION
The development of personalized neurorehabilitation and augmentation technologies implies the profound understanding of the mechanisms underlying an individual's motor ability and impairment.The development of 1) clinically viable techniques for interfacing with the human nervous system in vivo and 2) real-time subject-specific neural data-driven musculoskeletal models, are key factors for achieving this understanding.The reliable determination of musculoskeletal forces as a function of an individual's neural drive to muscles will provide the foundation to understand the neuromechanical interplay underlying in vivo movement function, pathology, and recovery.This will facilitate the design of personalized surgical and neurorehabilitation interventions, will inform the design of biologically inspired limbs, and will enable establishing neuromuscular HMIs for replacing or extending an individual's neuromusculoskeletal function.

Fig. 4 .
Fig. 4. Schematic structure of the open-loop neural data-driven modeling formulation.Estimates of neural excitations are derived from raw EMG data typically as linear envelopes, low-dimensional synergy modules, or motor unit spike trains [see Section III-A, Fig. 3(B)].The neural activation component (A) converts input neural excitations into neural activation using a second order muscle twitch model and a nonlinear transfer function [65], [98].The MTU kinematics component (B) synthetizes the MTU paths defined in the subject-specific geometry model (see Section III-C) into a set of MTU-specific multidimensional cubic B-splines.Each B-spline computes MTU length and moment arms as a function of input joint angles[125].The MTU dynamics component (C) solves for the dynamic equilibrium between muscle fibers and series elastic tendons in the production of MTU force[7],[176].It employs a Hill-type muscle model informed by MTU length and neural activations from the previous two components.The joint dynamics component (D) transfers MTU forces to the skeletal joint level using MTU moment arms.In the offline calibration component (E) initial nominal parameters are repeatedly refined, as part of a least-squares optimization procedure, so that the mismatch between model's predicted joint function and the experimentally recorded joint function is minimized[7],[8],[65].