Multibody dynamics-based musculoskeletal modeling for gait analysis: a systematic review

Abstract

Beyond qualitative assessment, gait analysis involves the quantitative evaluation of various parameters such as joint kinematics, spatiotemporal metrics, external forces, and muscle activation patterns and forces. Utilizing multibody dynamics-based musculoskeletal (MSK) modeling provides a time and cost-effective non-invasive tool for the prediction of internal joint and muscle forces. Recent advancements in the development of biofidelic MSK models have facilitated their integration into clinical decision-making processes, including quantitative diagnostics, functional assessment of prosthesis and implants, and devising data-driven gait rehabilitation protocols. Through an extensive search and meta-analysis of over 116 studies, this PRISMA-based systematic review provides a comprehensive overview of different existing multibody MSK modeling platforms, including generic templates, methods for personalization to individual subjects, and the solutions used to address statically indeterminate problems. Additionally, it summarizes post-processing techniques and the practical applications of MSK modeling tools. In the field of biomechanics, MSK modeling provides an indispensable tool for simulating and understanding human movement dynamics. However, limitations which remain elusive include the absence of MSK modeling templates based on female anatomy underscores the need for further advancements in this area.

Introduction

Gait or human walking signatures have long been used to quantify underlying health conditions, ranging from neurological and musculoskeletal (MSK) conditions to cardiovascular and metabolic disease, and to ageing associated ambulatory dysfunction and trauma [1]. Conventionally, gait is analyzed using measured anthropometric, kinematic, external kinetic, and muscle activity parameters. With the advancement in computing power, solving complex equations has been made feasible leading to the development of multibody dynamic models and systems mimicking humans [2].

Several multibody systems have been used in MSK modeling to understand and simulate complex interactions within the human body and of the human body with the environment. These include Empirically-based Multibody Dynamics [3], fast ligament models [4], OpenSim [5], Combined Multibody Musculoskeletal Dynamic Modeling and Finite Element Modeling [6] and commercial multibody modeling MSK systems like the AnyBody Modeling System (AMS) [7], Human Body Model [8], Software for Interactive Musculoskeletal Modelling [9], and Biomechanics of Bodies [10]. These software packages play a crucial role in various applications, including biomechanical evaluation of healthy and impaired individuals.

Gait analysis coupled with MSK modeling emerges as a state-of-the-art approach to delve deeper into the biomechanical intricacies of gait. This comprehensive method integrates multimodal data from motion capture systems for kinematics; force plates or instrumented walking mats/treadmills for ground reaction forces and moments (GRFs/GRMs); and muscle activity from electromyography (EMG) systems. The multidimensional data is then fed into MSK modeling software, which uses a generic MSK model to calculate the intersegmental joint reaction forces (JRFs), joint reaction moments (JRMs), as well as muscle forces using either inverse or forward dynamic analysis.

The primary aim of this review is to provide a systematic literature overview of different existing multibody MSK modeling methodologies, and their applications for gait analysis. More specifically, the study summarizes experimental protocols and data required to run MSK simulations; describe the most commonly used MSK modeling templates; highlight the personalization of these generic templates to match subject anthropometrics; outline the commonly employed solvers i.e., inverse kinematics (IK), inverse dynamics (ID), and identify the muscle recruitment techniques to compute joint kinematics, kinetics, muscle forces, and external kinetics. Moreover, this article also describes post processing of gait data for event detection, and the importance of joint JRFs, JRMs, and muscle parameters, as well as relevant applications of MSK modeling. This review is organized as follows: Sect. 2 introduces the adopted methodology, including the search criteria and meta-analysis. Section 3 reviews the conventional parameters quantifying gait and significance of MSK modeling in gait analysis. Section 4 describes the most common kinematic, kinetic, and muscle activity measurement techniques used for model input data, while Sect. 5 elaborates on current state-of-the-art platforms used in MSK modeling. Sections 6 and 7 describe generic MSK modeling templates and their personalization techniques, respectively. Section 8 discusses solving formulations, including inverse kinematics, dynamics, and muscle recruitment. Section 9 highlights post processing of gait data, while Sect. 10 summarizes relevant applications.

Review methodology

The process of finding studies, evaluating them for suitability, and extracting data for this review meticulously adhered to the guidelines outlined in the Preferred Reporting Items for Systematic Reviews and Meta-Analysis (PRISMA) statement [11] as presented in Fig. 1.

Search criteria

A keyword search was performed in PubMed, Scopus and Web of Science databases using a combination of search terms, including Musculoskeletal Modelling, Gait Analysis, Gait Joint Kinematics, Gait Joint Kinetics, OpenSim, and AnyBody. Filters applied to these platforms included full text, clinical trial, comparative study, technical report, validation study, and English language. The reference lists of the included articles were further reviewed to identify any relevant publications that met the predefined inclusion criteria.

Fig. 1

Diagram illustrating PRISMA article selection process

Inclusion and exclusion criteria

The inclusion and exclusion criteria for this study were established based on its objectives. The primary goals included describing the methodology for performing musculoskeletal (MSK) modeling in gait analysis, highlighting key equipment/manufacturers and software used, identifying commonly employed MSK models, describing methods for subject-specific personalization, summarizing approaches for solving statically indeterminate problems, detailing post-processing techniques, and showcasing practical applications of MSK modeling tools. According to these objectives, this systematic review includes both original research articles and review papers which met the following inclusion criterion: (i) studies centered on instrumented gait analysis employing MSK modelling; (ii) research assessing the validity and reliability of gait analysis; (iii) studies introducing MSK models; (iv) comparative studies of various MSK models; (v) research involving the personalization of MSK models; and (vi) investigations showcasing diverse applications of gait analysis. Following a careful screening process, 116 articles matched the established inclusion criterion and were subsequently included in the review.

Categorization and meta-analysis of studies

The selected studies were further primarily classified into MSK modeling approaches and applications. Among 116 selected publications, eleven studies focused on motion capture and GRF prediction; ten centered on MSK model scaling and personalization; fifteen compared MSK models or software; fourteen introduced custom models and plugins; and fourteen described gait analysis outcomes. The remaining fifty-two research studies explored the use of gait analysis to investigate various disorders and assess the efficacy of pertinent surgical procedures.

A meta-analysis of the selected studies was also performed, and the results are presented in Table 1. The meta-analysis conducted indicates that optical motion capture (MoCap) systems, typically comprising 6 to 10 cameras, are predominantly used to measure kinematics, with Vicon being the leading manufacturer of these systems. Additionally, most studies employed force plates operating at frequencies between 100 and 2000 Hz to capture GRFs, with AMTI, Bertec, and Kistler as the primary manufacturers. This meta-analysis also identifies the most commonly used software programs, MSK models, scaling methods, and the types of impairments analyzed using MSK modeling. The relevant statistics from the meta-analysis are highlighted in the respective sections of the study.

Table 1 Meta-analysis of the studies (FP – force plate; TM – treadmill)

Parameters quantifying gait and significance of MSK modeling

Changes in the signature of gait, as shown in Fig. 2, or the unique sequential walking pattern in humans, reveal key information about the status and progression of numerous underlying health challenges. These range from neurological and MSK conditions, to cardiovascular and metabolic disease, and to ageing associated ambulatory dysfunction and trauma [1]. Gait assessment typically comprises of the quantification of spatiotemporal characteristics (cadence, cycle time, speed, stride length, etc.), joint kinematics (joint angular positions/velocities/accelerations), external and internal joint kinetics (joint reaction forces and moments, ground reaction forces and moments), and muscle activity. Further details on parameters quantifying gait can be found in [1, 12].

Fig. 2

Difference phases of gait

Conventionally measured gait parameters, including spatiotemporal parameters, joint kinematics, center of mass trajectories, and GRFs and GRMs, are insufficient to quantify intersegmental joint loading, which is critical for biomechanical assessment. Joint moments, in essence, signify the tension exerted on structures like ligaments or muscles surrounding a joint. These moments are produced when the bones restrict the joint’s range of motion (ROM). A certain amount of power is absorbed or stored during ligament stretching, and this power is subsequently released during power generation [13].

Inverse dynamic analysis, using measured GRFs, computes external joint kinetics but overlooks the contribution of load-sharing components, such as muscles, ligaments, and the joint structure itself. Multiple pathologies including joint deformities, injury, sensory deficiencies, or impaired motor control, disrupt the coordinated patterns between muscle and joints. In such cases, accurately determining the damage to the affected muscles and ligaments is critical, which requires a comprehensive understanding of muscle and internal joint loading, in addition to conventionally measured gait parameters [14].

Utilizing MSK modeling coupled with inverse dynamics provides a powerful tool to comprehensively evaluate and analyze gait patterns, including synergies between muscle activation, joint kinetics, and motion. These advanced models allow medical professionals to quantitatively assess deviations from normal gait, facilitating the development of data-driven individualized treatment strategies and establishing a framework for monitoring patients’ progress over time [14].

MSK models also play a very important role in advancing prosthetics, exoskeletons, and computing loading for implants. This simulation-driven approach encourages the design of devices and solutions which mimic natural biomechanics, enhancing comfort, functionality, and overall user satisfaction. Moreover, MSK models facilitate the development of specialized finite element simulations, especially for complex joints such as the knee. The subsequent sections detail the methodology for conducting gait analysis using MSK modeling.

Kinematic, kinetic, and muscle activity data capture for MSK modeling

The primary sources of input for gait analysis using MSK modeling comprise of kinematic motion capture data, external kinetics represented by GRFs, as well as muscle activation data for some models.

Kinematic motion capture

The field of motion capture has undergone significant advancements in recent decades, introducing various techniques such as optical motion capture, wearable sensors, and fluoroscopic method to accurately measure human kinematics [15].

Optical motion capture (MoCap) systems are broadly classified into marker or marker-less systems. Marker-based systems further branch into active (using Light Emitting Diodes (LEDs) to indicate human body segments for MoCap) and passive marker-based (employing retroreflective markers applied to human skin) MoCap systems. Various protocols, including Plug-in-Gait, Modified Plug-in-Gait, Helen Hayes, Modified Helen Hayes, CAST model, and Cleveland Clinic marker sets, dictate the attachment of markers to subjects for motion assessment.

However, the application of markers is limited by potential shifts in underlying bones, leading to errors in joint kinematic estimation, known as skin tissue artefact (STA). While the kinematics obtained from these systems can serve as input when general output patterns are desired, caution is advised when more precise estimates are required [16].

The meta-analysis indicates that 86.21% of studies utilized optical MoCap systems for collecting kinematic data, with the number of cameras employed ranging between 4 and 16, predominately using 8 (22%), 10 (23%), or 12 (10%) cameras. Camera frequencies varied between 30 and 250 Hz, with the majority using frequencies of 100 Hz (40%), 120 Hz (9%), and 200 Hz (10%).

Well-known producers of optical MoCap systems include Vicon [17], Qualisys [18], Motion Analysis Corporation [19], and BTS Bioengineering [20]. Optical systems, despite their common usage in MoCap, come with significant drawbacks including high cost, susceptibility to STA, requirement of special expertise/training for operation, and system set up, tedious marker placement procedures, and calibration challenges. These limitations have restricted the widespread application of optical MoCap systems in clinical settings and inspired the development of marker-less alternatives.

Motion tracking in marker-less MoCap systems typically relies on three-dimensional depth cameras. Microsoft pioneered these cameras with the introduction of the Kinect in 2010, and subsequently, two other versions: Kinect v1 and Kinect v2. The Kinect v1 camera employs a structured light sensor with a resolution of 320 × 240 pixels, generating a phase-shifted infrared pattern to capture depth information through distortion. However, it has several drawbacks, including the need for dual viewpoints for proper depth capture, sensitivity to ambient light, and depth discontinuities. In contrast, the Kinect v2 overcomes these limitations by utilizing a 0.13 μm system-on-chip time-of-flight sensor with a spatial resolution of 512 × 424 pixels. Depth information is derived from the phase difference between emitted and received signals, achieved by emitting infrared light onto the scene. Unlike its predecessor, time-of-flight sensor-based cameras require only one perspective for depth data acquisition and exhibit resilience to ambient light. Despite being discontinued since 2017 [15], Kinect depth cameras continue to find applications in motion analysis studies [21,22,23].

The wearable sensor used for motion capture, as shown in Fig. 3, uses an Inertial Measurement Unit (IMU) to measure angular motion and linear acceleration of a body segment in three dimensions, employing a gyroscope and an accelerometer, respectively. The kinematic data form IMUs can be used by MSK modeling software like AMS or OpenSim [24,25,26,27] to predict joint reaction and muscle forces. Bailey et al. [25] conducted an assessment of the IMU-driven MSK model’s sensitivity and validity, concluding its suitability for evaluating joint angle time series, variability magnitude, ROM, and dynamic stability. Additionally, Karatsidis et al. [26] investigated the validity of measuring kinematic data using IMUs. Xsens stands out as a prominent manufacturer of IMU sensors for MoCap [28].

Fig. 3

Motion capture using inertial measuring units and instrumented treadmill during gait

Another promising technique for measuring joint kinematics without the limitations of soft tissue artifact (STA) is fluoroscopy. Fluoroscopy is a medical imaging method used to capture continuous high-speed X-ray images of the human body. The technique typically uses a C-shaped arm to align the X-ray source and detector. Recent advancements in fluoroscopy have enabled the monitoring of larger fields of view (up to 0.5 m × 0.5 m [29]) at high frequencies (typically up to 100 Hz [30]) while maintaining a radiation dose equivalent to less than a month’s worth of background radiation [30].

Fluoroscopic kinematic measurement involves capturing two-dimensional images of bones, which are then used to determine three-dimensional orientations by registering a separately obtained personalized 3D model to the images. This approach allows for the accurate measurement of six degrees of freedom (three translational and three rotational movements) of the joints. Fluoroscopy is often combined with instrumented treadmills for analysing walking and running [31].

There are two types of X-ray fluoroscopic techniques: mono-planar (using a single fluoroscope) and bi-planar (using two fluoroscopes). Mono-planar fluoroscopy has been shown to have an average root mean squared error (RMSE) of 0.53 mm for in-plane translations, 1.59 mm for out-of-plane translations, and 0.54° for rotations [32]. Bi-planar fluoroscopy improves out-of-plane translation accuracy, with an RMSE of 0.69 mm [33].

While fluoroscopy is widely used [34,35,36] for the assessment of joint kinematics, it has several limitations. First, the technique exposes subjects to radiation, necessitating its use by trained personnel in radiation-controlled environments. Additionally, the field of view is limited to a few joints, unlike conventional optical MoCap systems or IMUs, which can measure full-body kinematics over a larger volume. Finally, the cost and time required to process fluoroscopic images currently limit its feasibility for widespread use in clinical settings [31].External Kinetics Measurement

Ground reaction forces, moments, and center of pressure are typically measured using force plates and instrumented treadmills. However, force plates present several limitations. Firstly, subjects may need to alter their natural gait to ensure contact with the force plate’s limited area. Secondly, capturing both feet and multiple steps often necessitates the use of multiple force plates, which is costly. Thirdly, incorporating force plates in motion analysis can lead to dynamic inconsistencies, resulting in residual forces and moments during ID computations [27].

Alternatively, instrumented treadmills equipped with force measurement sensors can also be employed to measure GRFs. However, instrumented treadmills share similar drawbacks to force plates, including dynamic inconsistency and the potential alteration of the subject’s natural gait pattern. Among the 116 studies reviewed, 73 (62.93%) utilized force plates, while 25 (21.55%) employed treadmills to measure GRFs. Notable manufacturers of force plates and treadmills include AMTI [37], Bertec [38], and Kistler [39]. Both force plates and treadmills typically operate at a frequency of 1000 Hz.

It is worth noting that AMS employs an algorithm to predict GRFs along with muscle forces, as explained by Fluit et al. [40]. Although initially tailored for optical MoCap systems, recent adaptations have extended its applicability to depth cameras [21, 22] and IMUs [26]. The overview of studies centered on motion capture and GRF predictions is detailed in Supplementary Table S1.

Muscle activity measurement

Muscle activity is typically evaluated using wireless surface EMG electrodes. Leading manufacturers of EMG sensors include Delsys [41], Noraxon [42], and BTS Bioengineering [43]. Across the 116 selected studies, the number of EMG sensors utilized ranged from 1 to 16, contingent upon the specific muscles being studied. However, the use of 8 and 16 sensors predominated for collecting EMG data. Additionally, most EMG sensors operate at a frequency of 1000 Hz to capture muscle activity accurately. This muscle activity data is either integrated with kinematic and GRF data as input to the MSK model for predicting muscle forces or employed to validate the model by comparing muscle activations.

Computational musculoskeletal modeling tools

The meta-analysis of the selected studies indicated that OpenSim was used in 82.76% of the studies, while AMS was employed in 14.66%. Hence, this section exclusively describes these software packages.

OpenSim

OpenSim, an open-source software freely available to users, is specifically designed for motion analysis, modeling, and simulation of MSK systems [5]. It offers robust capabilities for numerical integration and solving constrained non-linear optimization problems. The procedure for conducting MSK simulations in OpenSim is presented in Fig. 4 and includes the following steps: 1- Capturing motion data using a motion capture system (as explained in section 4); 2- Selecting an appropriate MSK model (described in section 6); 3- Scaling the model to match the subject’s anthropometric dimensions (refer to Sect. 7); 4- Executing IK analysis; 5- Conducting ID analysis (further elaborated in Sect. 8); 6- Running Residual Reduction Algorithm (RRA); 7- Computing muscle forces by solving an optimization problem; and 8- Post-processing to derive the desired outputs (detailed in section 9).

Fig. 4

Flowchart illustrating the process of performing an MSK simulation using OpenSim

OpenSim is versatile in accommodating motion capture data from various sources, including MoCap, and IMUs. Moreover, it offers a wide array of muscle recruitment strategies, including static optimization, computed muscle control (CMC), EMG-driven, and EMG-informed optimizations. The RRA in OpenSim is designed to tackle dynamic inaccuracies arising from noise, MoCap errors, and model assumptions. These inaccuracies can lead to deviations from Newton’s second law i.e., $\:\varvec{F}\ne\:m\varvec{a}$, which the RRA aims to rectify by adjusting the mass and center of mass of body segments [44]. One of the distinctive advantages of OpenSim is its flexibility, allowing users to customize and generate their own models to suit specific research needs and objectives.

AnyBody musculoskeletal modeling system

AMS is a commercial closed-source software developed by AnyBody Technology [7] for MSK modeling and simulation. Figure 5 depicts a flowchart illustrating the procedural steps involved in conducting an MSK simulation through AMS.

While the simulation steps in AMS resemble those in OpenSim, the primary distinction lies in AMS’s omission of RRA processing; instead, it proceeds directly to muscle computation and post-processing. AMS also accommodates kinematic data from diverse motion capture technologies, including MoCap, IMUs and depth cameras. One notable advantage of AMS is its capability to predict GRFs during simulation, thereby eliminating the need for them as input data. The MSK models accessible within AMS include the Twente Lower Extremity Model (TLEM), TLEM 2.0, and GaitFullBody.

Fig. 5

Flowchart illustrating the process of performing an MSK simulation using the AnyBody Modeling System

Comparison between AnyBody and OpenSim

Table 2 provides a qualitative comparison between AMS and OpenSim software packages. While this table addresses qualitative aspects, quantitative comparisons have been conducted by Alexander et al. [45] and Trinler et al. [46]. Alexander et al. [45] investigated the influence of walking velocity on muscle force and joint loading estimations using both OpenSim and AMS. Their findings indicated a high correlation coefficient (> 0.96) between joint loading predictions made by both software packages across all speed scenarios. Additionally, muscle force correlation coefficients exceeded 0.7 in 80% of speed conditions, except for the tibialis posterior and biceps femoris long head muscles.

Trinler et al. [46] conducted a comparison of muscle forces, kinetics, and kinematics estimations between AMS and OpenSim during gait. Notably, discrepancies were consistently observed in ankle and hip kinematics in the sagittal plane throughout the gait cycle. Specifically, the ankle joint exhibited the most significant difference between OpenSim and AMS predictions during 50% of the gait cycle. These disparities were attributed to variations in anatomical and anthropometric definitions between the two software platforms. Consequently, it is imperative to consider these differences when conducting and interpreting gait analyses.

Table 2 Qualitative comparison of AMS and OpenSim software

Generic musculoskeletal model templates

Derived from cadaveric data, generic MSK models primarily vary in degrees of freedom (DOF), joint kinematics, bone geometries, coordinate system definitions, and muscle-tendon unit parameters [47]. Over recent decades, several generic MSK models with varying levels of complexity have been developed and documented in the literature. The following section provides a review of the generic MSK models available in OpenSim and AMS.

Musculoskeletal models in OpenSim

In OpenSim, lower limb models include the Delp model [48], Gait2392, Gait2354, Arnold Two Leg Model (ATLM) [49], and London Lower Limb Model (LLLM) [50]. Meanwhile, full body models include the Full Body Running Model (FBRM) by Hamner et al. [51], Rajagopal model [52], Lai Model [53], and Full Body Model 2016 (FBM16) by Caruthers et al. [54]. The findings from the meta-analysis reveal that the Gait2392 and Gait2354 models were utilized in 42.24% of the studies, while 12.07% employed the Delp Model and ATLM. Additionally, the Rajagopal model was used in 9.48% of the studies. Other models such as LLLM, FBRM, Lai Model, FBM16, etc. collectively contributed to approximately 8.62% of the selected research studies.

The lower limb MSK model, developed by Delp et al. [48], was designed to investigate the impact of surgical modifications to bone geometries and muscle parameters on the function of muscles. This generic model has 43 muscle-tendon actuators, 7 DOF, 7 body segments (pelvis, femur, patella, tibia, talus, foot, and toes), and represents a 180 cm tall male. The muscle-tendon actuator units in this model are constructed by scaling a generic Hill-type model using muscle parameters obtained from studies conducted by Wickiewicz et al. [55] and Friederich and Brand [56]. Physiological cross-sectional area (PCSA) values for muscles, as published by Friederich and Brand [56] and Wickiewicz et al. [55], respectively, were obtained through experiments involving young and elderly cadavers. To ensure data correlation with young adults, the values derived from elderly cadavers were appropriately scaled within the model by Delp et al. [48].

The Delp model [48] underwent enhancements, evolving into the Gait2392 model, aimed at predicting internal loads within MSK systems [57, 58], and investigating individual muscle function during movement [59]. This model represents a male weighing 75.16 kg and standing 180 cm tall in its unscaled function, and has 92 muscle-tendon actuator units, 12 segments, and 23 DOF in the lower extremities and torso. The Gait2392 model incorporates the planar knee model from Yamaguchi and Zajac [60], the anthropometric and low back joint from Anderson and Pandy [61], as well as the lower body joint definition from Delp et al. [48]. Apart from the toes and hindfeet, the inertial and mass parameters of ten of the twelve bone geometries or body segments are derived from Anderson and Pandy [61]. On the other hand, the inertial and mass parameters for the toes and hindfeet are established by defining a volume created by joining the vertices or surface coordinates of a size-10 tennis shoe. The mass of each segment is then calculated by numerically integrating density over the foot volume, assuming a uniform density of 1.1 gcm^− 3. Muscle parameters are obtained from Wickiewicz et al. [55] and Friederich and Brand [56], similar to the Delp et al. [48] model. However, scaling factors are updated to increase muscle strength when converting muscle data from the elderly to young adults. Additionally, Gait2354, a simplified version of the Gait2392 model, has 54 muscle-tendon actuator units and 23 DOF. The primary objective behind reducing the number of muscles is to enhance computational speed, particularly for learning and demonstration purposes.

Hamner et al. [51] developed the Full Body Running Model (FBRM) to investigate the role of muscles in propelling the body’s center of mass forward and upward during running. Featuring 92 muscle-tendon actuator units, 12 rigid body segments, and 29 DOF, the FBRM represents an enhanced and refined iteration of the Gait2392. To transform the Gait2392 into the FBRM, the arms were integrated as idealized torque actuators. Similar to the Delp et al. [48] and Gait 2392 models, muscle parameters for the FBRM are sourced from Wickiewicz et al. [55] and Friederich and Brand [56].

Arnold et al. [49] introduced the Lower Limb Model 2010, also known as ATLM, as an extension of the Delp et al. [48] and Gait2392 models, aimed at predicting muscle forces and fiber lengths in the lower body. This model uses 44 Hill-type muscle-tendon units as actuators in the lower extremity, with a total of 23 DOF. The default version of the ATLM in OpenSim contains 14 body segments and represents an adult male height of 170 cm in height and 75 kg in weight. Derived from Arnold et al. [62], the rigid body segments of the ATLM underwent significant improvements over the Delp et al. [48] model. These enhancements include revisions to the pelvis coordinate system, re-digitalization of the femur, re-orientation and rescaling of the fibula, and replacement of the tibia. Mass and inertial parameters for the body segments are sourced from Anderson and Pandy [61], while properties for the hindfoot and toes are computed by numerically by integrating a constant density of 1.1 gcm^− 3 over the segment volume. The ATLM incorporates data from Ward et al. [63], gathered from 21 cadavers (12 female and 9 male) with an average age of 83 ± 9 years, to include muscle parameters. This extensive dataset enables the ATLM to be considered as a generic lower extremity model. The model primarily defines muscle-tendon unit pathways using ellipsoidal wrapping surfaces, thereby enhancing its ability to accurately predict fiber lengths and muscle forces. However, this approach increases dynamic simulation time compared to simple models, such as Gait2392 and Gait2354. Despite its widespread use in simulating gait, the ATLM has two notable drawbacks in muscle-driven dynamic simulations [16, 47, 64,65,66,67,68,69]. Firstly, its extensive use of ellipsoidal wrapping surfaces makes simulations computationally expensive, rendering it impractical for forward muscle-driven simulations. Furthermore, the incorporation of muscle data from older cadavers with significant muscular atrophy results in inaccurate muscle force predictions for young healthy adults.

To address these challenges, Rajagopal et al. [52] created a full body model towards facilitating MSK modeling of human gait. As illustrated in Fig. 6(a), the Rajagopal model has 22 rigid body segments and 37 DOF, with 20 in the lower extremities and 17 in the upper extremity. The unscaled version of the model represents a male individual weighing 75 kg and standing 170 cm tall. Incorporating enhancements over previous models, the Rajagopal model comprises of 17 torque actuators in the upper extremity and 80 Hill-type muscle-tendon units as actuators in the lower extremity, as depicted in Fig. 6(b). The top extremity of the model is derived from the Hamner et al. [51] model, while the lower extremity is generated from ATLM. Noteworthy augmentations include updated muscle wrapping surfaces and peak isometric muscle forces as compared to ATLM and prior models. Muscle parameters, such as ideal fiber length, pennation angle, and tendon slack length, are obtained from the cadaveric investigation by Ward et al. [63]. On the other hand, the peak isometric forces of muscles are determined using muscle volume data from Handsfield et al. [70]. This dataset involves magnetic resonance imaging (MRI) measurements of muscle volume from 24 young healthy subjects, including 8 females and 16 males, with an average age, height, and weight of 25.5 ± 11.1 years, 171 ± 10 cm, and 71.8 ± 14.6 kg, respectively. To streamline the computations, the Rajagopal model modifies the elliptical wrapping surfaces of muscle paths from ATLM to cylindrical wrapping surfaces, leveraging an existing analytical solution for muscle wrapping around cylindrical surfaces. This adjustment minimizes computational costs. An example of muscle wrapping around a cylindrical surface is illustrated in Fig. 6(d), while Fig. 6(e) depicts the modeling of the patellar ligament in the Rajagopal model.

Fig. 6

Rajagopal full body MSK model: (a) an illustration displaying 23 DOF and 22 rigid bodies, (b) a representation of 80 Hill-type muscle-tendon actuator units that are integrated into the lower body structure, (c) application of many distinct muscles to replace those muscles with large attachment areas, (d) wrapping of muscle surface around cylindrical surfaces, and (e) modeling of patellar ligament [52]

Lai et al. [53] refined the Rajagopal model to simulate running and pedaling in addition to walking. Fundamental attributes, such as the number of rigid body segments, DOF, muscle-tendon units, and sources of muscle parameter, remain consistent with the Rajagopal model. However, refinements are introduced to the muscle-tendon parameters and pathways of 22 muscles without altering the maximum isometric muscle forces. In the Lai model, adjustments are made to knee joint kinematics to simulate actions involving high knee flexion more accurately, such as pedaling. Notably, the range of motion of the knee is extended from 0-120º to 0-140º, and the origin-to-insertion pathways of knee muscles are updated. Another important refinement involves the motion of the tibia relative to the femur. Utilizing mathematical relations similar to those in the ATLM and Rajagopal models, proximo-distal and anterior-posterior translations, abduction-adduction, and internal-external rotations are determined as functions of the knee flexion-extension angle. Of particular note, the Lai model includes the medio-lateral translation of the tibia relative to the femur, enhancing its fidelity in capturing complex knee movements during various activities.

Caruthers et al. [54] developed the FBM16, a full body MSK model designed to investigate muscle forces and their associated accelerations during sit-to-stand exercises in healthy young adults, with applications in gait simulations [67]. The FBM16 integrates components from four distinct models: the lower extremity model is based on the ATLM by Arnold et al. [49]; the upper extremity is based on the lumbar spine model by Christophy et al. [71]; the upper limb model is based on Saul et al. [72]; and the head and neck MSK biomechanics model is based on Vasavada et al. [73]. Comprising 194 muscle-tendon units, 24 rigid body segments, and 23 DOF, the FBM16 offers a comprehensive framework for investigating complex MSK interactions. Given the inherent anthropometric discrepancies among the constituent models, the ATLM serves as the baseline, with the bodies of the other models scaled to match a male subject standing 170 cm tall and weighing 75.337 kg to ensure consistency across the integrated components. The masses, inertial characteristics, muscle-tendon parameters, and wrapping surfaces are also appropriately scaled using conversion ratios, as detailed in Caruthers et al. [54].

Modenese et al. [50] developed the LLLM, a unilateral MSK model of the lower extremity, to accurately predict hip contact forces, Unlike previous models combining data from disparate sources, the LLLM aims for coherence by using a single dataset for the entire MSK system. The LLLM is based on the anatomical dataset from Klein Horsman et al. [74], obtained from the right leg of a single male cadaver (aged 77 years, height 174 cm, weight 105 kg). This dataset provides comprehensive information crucial for MSK modeling, comprising of joint kinematics, geometry, muscle contraction characteristics, and precise muscle attachment locations [50]. The unilateral LLLM offers a focused yet comprehensive representation of lower limb dynamics, representing the lower limbs with 6 rigid bodies, 5 joints, and 163 muscle-tendon actuators. The 11 DOF inherent to this model provide realistic simulations capturing the intricacies of human movement. To enhance this model’s accuracy, Modenese et al. [50] added several adjustments to various aspects, including the talocrural joint axis movement, muscle pathways, and the insertion of distal bundles, based on insights collected from the dataset published by Klein Horsman et al. [74].

Musculoskeletal models in AnyBody modeling system

In the AMS, two lower body models stand out: the TLEM [74], as well as its updated version, TLEM 2.0 [75]. Additionally, AMS offers a comprehensive full body model known as GaitFullBody. Through the meta-analysis, it was found that 12.07% of the studies included in this review relied on either TLEM, TLEM 2.0, or GaitFullBody models for gait analysis.

Prior to the introduction of the TLEM, MSK models were typically developed using different input data sourced from various studies as appropriate. For example, anatomical geometry, muscle-tendon characteristics, attachment locations, and pathways were often integrated from multiple references, such as Friederich and Brand [56], Yamaguchi and Zajac [60], and Wickiewicz et al. [55]. Notably, no single anatomical investigation had comprehensively reported all necessary parameters until Klein Horsman et al. [74]. Klein Horsman et al. [74] work marked the inception of TLEM, a model meticulously developed from the detailed analysis of a single male specimen’s right lower extremity. The specimen, 77 years old, 105 kg in weight, and 174 cm in height, provided the foundational data for TLEM’s anatomical geometry, inertial parameters, joint kinematics, muscle contraction characteristics, and attachments locations [74]. Implemented within the AMS, TLEM has 12 rigid body segments, including the pelvis, HAT (Head, Arms, and Trunk), bilateral femurs, tibia-fibula, patella, foot, and talus. Its articulations include 11 joints: L5-S1, bilateral subtalar, talocrural, femur/patella, knee, and hip. Featuring 159 three-element Hill-type muscle-tendon units as actuators, the TLEM offers a total of 21 DOF [76, 77].

Notably, the accuracy of personalized MSK model predictions hinges on the precision of the scaling technique used to tailor the generic model to individual subjects. Among the most accurate methods is the medical image scaling technique, presented by Pellikaan et al. [78]. This approach involves morphing computerized tomography (CT) scanned images of the subject onto a pre-existing template with predefined muscles’ lines-of-action and attachment sites. Harnessing the power of medical image scaling, TLEM 2.0 is designed to deliver personalized MSK models with integrated templates for precise geometry. This entails precise calculations of the muscle-tendon line-of-action, moment arm, and overall length. Unlike its predecessor, TLEM 2.0 prioritizes personalized models over generic ones. The model’s development involved dissecting a male specimen, estimated to be 85 years old and weighing 45 kg. Prior to dissection, MRI and CT scans were performed on both lower limbs to acquire the detailed anatomical data, including muscle and subcutaneous fat volumes, as depicted in Fig. 7a. Further details on the measurement process can be found in Carbone et al. [75, 76]. Integrated within AMS, TLEM 2.0 comprises of twelve body segments, including the pelvis, HAT, biliteral femurs, tibia-fibula, patella, foot, and talus. The model’s articulations encompass eleven joints: knee, hip, femur/patella, talocrural, bilateral subtalar, in addition to L5-S1. Figure 7b illustrates TLEM 2.0’s configuration, featuring 21 DOF and 166 Hill-type muscle-tendon units as actuators [75, 76].

Alongside TLEM and TLEM 2.0, AMS introduces the holistic GaitFullBody model, a comprehensive representation of the MSK system. Its lower extremity inherits the meticulous design of TLEM by Klein Horsman [74], renowned from its anatomical fidelity and precise modeling parameters. Meanwhile, the lumbar spine segment draws from the insights of de Zee et al. [79], ensuring robust biomechanical accuracy in spinal dynamics. Similarly, the shoulder-arm complex is designed based on the research of the Delft Shoulder Group [80,81,82]. As highlighted by Skals et al. [83], GaitFullBody is considered as one of the most detailed generic full-body MSK models available, featuring a total of 42 DOF [84]. Supplementary Table S3 contains summary of the studies about the MSK models and OpenSim plugins.

Fig. 7

Twente Lower Extremity Model 2.0: (a) segmentation of the lower extremity including muscle and subcutaneous fat volumes, and (b) resulting MSK model with 12 rigid bodies, 11 joints, 166 muscle-tendon units, and 21 DOF in AMS [75, 76]

Comparison between generic MSK model templates

Several studies have investigated the impact of MSK model selection on the accuracy of predicting kinetics, kinematics, muscle activation patterns and forces. Roelker et al. [67] conducted a comprehensive analysis of four OpenSim models, namely Gait2392, ATLM, FBRM, and FBM16. Their findings revealed a correlation between model complexity and disparity between simulated and measured results, with Gait2392 demonstrating sufficient fidelity in capturing gait dynamics among young healthy adults. Lamberto et al. [16] extended this investigation, exploring the sensitivity of MSK models — ATLM, LLLM, and Gait2392 — to STA through kinematic and kinetic predictions using OpenSim. Despite variations in model parameters, the examined models exhibited consistent patterns in investigated variables, albeit with notable differences in magnitudes. In another study, Mathai and Gupta [66] evaluated the performance of four MSK models — Gait2392, LLLM, ATLM, and Hip2372 — in estimating hip joint forces in Total Hip Arthroplasty (THA) patients. Their analysis revealed minimum differences in error (less than 30% body weight (BW)) between measured hip JRFs and estimates from all models excluding ATLM. Additionally, during the late stance phase of gait, all models tended to overestimate hip JRF predictions, with LLLM demonstrating the most accurate predictions. Weinhandl and Bennett [69] assessed hip JRF predictions from four distinct OpenSim models (ATLM, LLLM, Gait2392, and hip2372) during gait. Their investigation identified LLLM as the optimal choice, exhibiting the highest correlation values and lowest Root Mean Square Error (RMSE) among the assessed models. Lastly, Curreli et al. [47] compared knee JRF predictions from three MSK models —ATLM, Rajagopal, and Lai — in a patient with Total Knee Replacement (TKR). This study concluded that while ATLM showed poor correlation, the Rajagopal model yielded the most accurate predictions compared to measured values, indicating its suitability for estimating knee JRFs. Supplementary Table S2 outlines the details of articles concentrating on software and model comparisons.

In summary, the selection of an MSK model should align with the specific objectives of the gait analysis. Gait2392 emerges as a suitable option for studies involving healthy young adults, whereas LLLM and Rajagopal models offer superior accuracy in estimating hip and knee JRFs, respectively. For a quick reference, a qualitative comparison of various MSK models is provided in Table 3.

Table 3 Summary of various MSK models used for gait analysis in OpenSim and AMS

Integrated MSK and finite element models

In addition to the MSK models in AMS and OpenSim, various research groups have developed custom workflows for computing muscle and joint reaction forces. One notable example is the concurrent finite element MSK framework developed by Wang et al. [85, 86]. In this framework, angular velocity and acceleration data from IMUs, or marker trajectories from optical MoCap systems, are input into a custom-built 22 DOF full-body MSK model to calculate joint angles through inverse kinematics (IK). Next, measured GRFs from force plates or instrumented treadmills and joint kinematics are used to compute joint moments, moment arms, and muscle lengths via Newton-Euler ID analysis. These computed variables are then fed into a custom MATLAB static optimization algorithm to predict muscle forces. The primary joint angles, GRFs, and muscle forces serve as boundary conditions for a knee finite element model, typically generated using CT scans and MRIs. The finite element model estimates secondary kinematics of the knee joint, which are fed back into the inverse dynamics (ID) analysis, static optimization, and knee contact mechanics calculations. This iterative process continues until the balance between joint moments and muscle forces converges [85, 86].

Marouane et al. [87] also integrated a detailed finite element model of the knee with an MSK model of the lower body to predict tibiofemoral contact forces and the location of contact centers on the lateral and medial plateaus. They reported that, during the early stance phase, the contact force was greater on the lateral plateau, while it was greater on the medial plateau during 25–100% of the gait cycle. Additionally, they computed larger excursions, exceeding 17 mm, in the location of the contact center. Similarly, Khatib et al. [88] incorporated a detailed finite element model of the knee into an MSK model to compute muscle and joint reaction forces in obese and normal-weight subjects during the stance phase of gait. They concluded that obese subjects experienced significantly higher forces compared to their normal-weight counterparts.

Furthermore, Xu et al. [89] investigated the impact of load carriage on the tibia during walking using an integrated MSK and finite element model. Full-body MoCap and GRF data were used to compute muscle and joint reaction forces under four different load-carrying conditions. These forces were then applied as loading conditions in a finite element model to calculate stresses and strains. They found that knee joint contact force increased by 26.2% when carrying a load equal to 30% of body weight compared to no-load conditions.

Halonen et al. [90] studied the effects of low wedge insoles (5° and 10°) and gait modifications (toe-in and toe-out) on the stresses produced in medial tibial cartilage. They used the output forces from an MSK model as input to a finite element model of the knee. Their findings revealed that the 10° insole reduced the second peak of knee load, while both toe-in and toe-out modifications reduced the first peak.

Shu et al. [91] proposed a unique subject-specific MSK model combined with a high-accuracy concurrent finite element knee model to study knee kinematics and mechanics. The maximum axial force acting on the knee joint was found to be 2.89% of body weight at 45% of the gait cycle. Additionally, they reported that approximately 65.7% of the axial knee joint force was borne by the medial cartilage.

In addition to knee models, several studies have combined finite element models of the hip [92,93,94,95,96,97] and foot [98,99,100] with MSK models. For example, Seo et al. [92] used joint reaction forces computed by an MSK model as input to a finite element model of the right femur to enhance the accuracy of the finite element analysis. Similarly, Ravera et al. [93] developed an integrated MSK and finite element model of the pelvis to study the mechanical behavior of tissues, finding good agreement between the hip joint contact forces predicted by both models.

Altai et al. [94] utilized muscle forces computed by an MSK model as boundary conditions in a finite element model of the femur to calculate femoral neck strains during gait, observing two peaks of maximum strain when the feet were in contact with the ground. Li [95] proposed a finite element MSK model of the lower body, incorporating the hip joint as a deformable contact model. This coupled model took joint kinematics and GRFs as input and performed an optimization problem, minimizing the sum of muscle forces squared. The model successfully computed muscle and contact forces, and stresses within the hip joint, showing good agreement with in vivo data.

Gaffney et al. [96] also used hip joint reaction forces as input to a finite element model to compute hip joint cartilage contact pressure during gait, finding reasonable agreement with predictions made by an elastic foundation contact model. Xiong et al. [97] used hip joint reaction forces computed by an MSK model as loading and boundary conditions in a finite element model of the hip joint, with maximum contact stress reaching 6.9 MPa in the acetabulum.

In the context of foot modeling, Akrami et al. [99] used ankle joint reaction forces, muscle forces, and GRFs as input to an MRI-based finite element model of the foot. A sensitivity analysis revealed that the model was significantly influenced by foot orientation. Kamal et al. [100] employed variables obtained from an MSK model as input to a finite element model of the foot, computing plantar pressure during gait and observing the highest concentration of 0.48 MPa in the forefoot region.

In summary, MSK and finite element models complement each other and can be effectively utilized to study soft tissue forces and improve implant design. However, these models have several limitations. Firstly, they are computationally intensive, requiring significant resources. Secondly, they are typically restricted to a single joint of interest, such as the hip, knee, or foot, limiting their broader applicability. Thirdly, this approach necessitates personalized finite element models created from MRI or CT scans, making the method costly and less accessible. Lastly, direct validation of these models remains limited in most cases, reducing their reliability for broader clinical use.

Neuromusculoskeletal models

MSK models available in tools, such as AMS and OpenSim, primarily rely on ID analysis to estimate muscle forces and activations based on observed human motion. However, these models are limited in their ability to identify or quantify the underlying causes of movement, such as the role of central pattern generators (CPGs) within the human neural system. To address this limitation, neuromusculoskeletal models that utilize forward dynamics to simulate human gait have been developed. Unlike ID-based models, these models generate movement through mathematical relationships between neural inputs and the human body’s response. In these models, CPGs activate muscles, which then drive the movement of body segments, effectively producing locomotion. To provide a comprehensive overview and showcase the latest advancements in MSK models, a representative list of neuromusculoskeletal models is presented in this section.

One of the pioneering neuromusculoskeletal models was presented by Taga [101], which included 20 muscles, 8 body segments, 7 pairs of neural oscillators for rhythm generation, and mechanisms for motor and sensory signal processing. This model, with over 100 parameters [102], successfully produced walking in a two-dimensional MSK model at speeds ranging from 0.7 m/s to 1.2 m/s by varying only a constant neural input to the oscillators. Hase et al. [103] later developed a more complex model comprising 60 muscles, 14 body segments, 16 pairs of neural oscillators, and 19 degrees of freedom (DOF). This model, with 125 parameters [102], used evolutionary computation to simulate both walking and running in a three-dimensional MSK model. Another notable contribution is the work of Günther and Ruder [104], who created a two-dimensional neuromusculoskeletal model based on the λ model of the equilibrium-point hypothesis [105]. This model, containing 20 parameters, was able to simulate walking by adjusting the trunk’s reference angle and stretch reflex feedback gain. Similarly, Song and Geyer [106,107,108] introduced a reflex-based neuromusculoskeletal model with 22 muscles and 7 body segments [107]. By integrating motor control into a two-dimensional MSK model with 30 parameters, they demonstrated that walking at speeds of 2.4 m/s to 4 m/s could be achieved by modifying just 9 of the parameters [106]. Their model [107] also featured 82 parameters that could be optimized to simulate both running and walking in a three-dimensional MSK model.

In contrast, Neptune et al. [109, 110] developed a motor control model based on the muscle synergy hypothesis to simulate walking. This model utilized measured EMG data to determine five activation pulses that control movement. By optimizing the duration, magnitude, and onset of these pulses, they achieved walking simulations in both two-dimensional (26 muscles, 13 body segments, and 13 DOF) [109] and three-dimensional (76 muscles, 12 body segments, and 23 DOF) [110] MSK models. Jo and Massaquoi [111] also contributed to this field by developing a two-dimensional neuromusculoskeletal model with 7 body segments, simulating walking using the muscle synergy hypothesis. Their model, which incorporated a linear combination of four activation pulses, was capable of simulating walking at speeds ranging from 0.6 m/s to 1.4 m/s by adjusting pulse magnitude, gait cycle duration, and the centre of mass reference. Finally, Aoi et al. [102, 112] developed a two-dimensional model with 18 muscles and 7 body segments, also based on the muscle synergy hypothesis. This model, comprising 69 parameters, was designed to simulate walking.

In summary, these neuromusculoskeletal models highlight the critical role of CPGs in simulating human gait, with ongoing research focused on identifying optimal parameters for accurately modeling walking dynamics. These models represent a significant advancement in understanding the neural and biomechanical interactions involved in human movement, offering potential applications in both research and clinical settings. Personalization of Generic MSK Template

Following the selection of a generic MSK model, the subsequent crucial step in the modeling process involves scaling or personalizing the model to match the anthropometric dimensions of the individual subject. The customizable properties of the generic MSK models include:

• Mass and Inertial properties of each model segment.

• Geometrical properties, encompassing muscle insertion and pathway sites, joint centers, centers of mass, etc.

• Characteristics of muscles and ligaments, such as ideal fiber length, pennation angle, slack length, etc.

Each of the aforementioned properties must be adjusted to generate a fully personalized MSK model. This task presents challenges, as many of these properties necessitate costly experimental testing and cannot be determined immediately. The subsequent sections provide a brief discussion of scaling key properties of generic MSK models.

Mass and inertial properties

This category includes the mass of each segment, its center of mass, and the mass moment of inertia about three axes of the segment. Body density is primarily influenced by the distribution of muscle, bone, and fat across the body which varies among individuals. In spite of this variation, MSK models typically use regression equations to estimate the mass of each segment based solely on total body mass [113]. Similarly, the center of mass is scaled from the generic model using the same methodology as other model points. Lastly, employing an idealized segment shape, often a cylinder, the mass moment of inertia is scaled according to segment masses and lengths [113].

Geometric scaling

Geometric scaling includes personalization of points (such as the center of mass, joint center, muscle origin, insertion, and via points, etc.), vectors (including joint axis, wrapping surface vectors, etc.), and surfaces (wrapping surfaces) from generic to subject-specific geometry. Scaling points entails mapping them using the function $\:\varvec{s}:{R}^{3}\to\:{R}^{3}$ from their source to their target positions. Since scaling surfaces and vectors might pose challenges, they are typically treated as points and scaled in a similar manner to points.

The easiest way to scale generic models is using linear scaling law [114]:

$$\:\varvec{s}=\varvec{S}\varvec{p}+\varvec{t}$$

(1)

In Eq. (1), $\:\varvec{s}$ denotes the scaled point in local segment’s fixed coordinate system, $\:\varvec{p}$ is the unscaled point, $\:\varvec{S}$ is a $\:3\times\:3$ scaling matrix, and $\:\varvec{t}$ is the translation vector. Despite its inclusion in the linear scaling law, the translation vector $\:\varvec{t}$ is not utilized by any of the scaling methods [113]. The diagonal values of the scaling matrix $\:\varvec{S}$ are the scaling ratios along the three directions of the segment coordinate system. Three possible definitions for $\:\varvec{S}$ are possible: uniform, non-uniform, and orthogonal.

If the matrix $\:\varvec{S}$ is defined as a diagonal matrix with all scaling factors equal, it is termed uniform scaling, as shown in Eq. (2), where k is the scale ratio.

$$\:\varvec{S}=k\varvec{I}$$

(2)

In non-uniform scaling, the scaling ratios in different directions vary. Equation (3) gives the matrix $\:\varvec{S}$ for non-uniform scaling. The scale ratios along the local x, y, and z axes are indicated by the symbols $\:{k}_{x}$, $\:{k}_{y}$, and $\:{k}_{z}$ respectively.

$$\:\varvec{S}=\left[\begin{array}{ccc}{k}_{x}&\:&\:\\\:&\:{k}_{y}&\:\\\:&\:&\:{k}_{z}\end{array}\right]$$

(3)

Finally, if the coordinate axes used for scaling do not align with the segment’s local coordinate axis, matrix $\:\varvec{S}$ is represented by Eq. (4), where $\:\varvec{A}$ is the orthogonal rotational matrix and $\:\varvec{S}\varvec{{\prime\:}}$ is the non-uniform scaling matrix as described in Eq. (3).

$$\:\varvec{S}=\varvec{A}{\varvec{S}}^{\varvec{{\prime\:}}}{\varvec{A}}^{\varvec{T}}$$

(4)

The next challenge lies in calculating the scaling ratios. Rasmussen et al. [114] proposed a technique known as the length-mass scaling law for calculating these ratios. In this approach, the ratios are calculated based on the mass and length of the segment. Equation (5) is used to calculate the scaling ratio along the longitudinal axis, denoted as $\:{k}_{y}$. The lengths of scaled and unscaled segments are indicated by the characters $\:{L}_{S}$ and $\:{L}_{T}$, respectively. Equation (6) is used to compute the scaling ratios in the cross-sectional directions, where $\:{k}_{m}=\:{m}_{S}/{m}_{T}$ represents the mass ratio of scaled ($\:{m}_{S}$) and unscaled ($\:{m}_{T}$) segments.

$$\:{k}_{y}={k}_{L}=\frac{{L}_{S}}{{L}_{T}}$$

(5)

$$\:{k}_{x}={k}_{z}=\sqrt{\frac{{k}_{m}}{{k}_{L}}}$$

(6)

Given the scaling background, techniques can be categorized into three types. Type-I involves scaling methods where body segment lengths are externally estimated from measured anthropometric data. Type-II encompasses methods utilizing static or dynamic motion capture experimental trials, including linear scaling, anatomical landmark scaling, and kinematic scaling. Lastly, Type-III relies on medical images to obtain more accurate segment lengths and joint centers. Table 4 provides a concise overview of each technique, along with their respective advantages and disadvantages. For more detailed descriptions, readers are referred to Andersen et al. [113]. The meta-analysis of the selected studies in this review revealed that the Type-II scaling method, based on motion data, was utilized by 52.69% of the studies, while Type-I, relying on anthropometric measurements, and Type-III, utilizing medical imaging, were employed by 26.72% and 8.62% of the studies, respectively. Supplementary Table S4 presents the details of articles related to the scaling of MSK models.

Table 4 Summary of geometrical scaling techniques

Muscle-tendon unit properties

In MSK modeling, muscles are typically represented using a Hill-type model, which offers a mechanical depiction of muscle-tendon dynamics [115, 116]. This model comprises of four key components: the contractile element, tendon, parallel elastic element, and serial elastic element. A comprehensive description of a Hill-type model requires specifying numerous parameters [113], but for personalizing a standard Hill muscle, four parameters are particularly significant. These include the optimal fiber length (OFL), maximal isometric force (MIF), pennation angle, and tendon slack length (TSL).

During the personalization process, OFL and TSL are typically linearly scaled alongside the geometric scaling of muscle attachment points to maintain consistent ratios (OFL/TSL). However, MIF is generally not subjected to scaling in MSK modeling [117]. Nevertheless, linearly scaling OFL and TSL while omitting MIF scaling can significantly affect the accuracy of muscle force predictions. A comprehensive review of various scaling techniques for Hill-type muscles can be found in Heinen et al. [118].

Inverse kinematics, inverse dynamics, and muscle recruitment

Following the personalization of the generic MSK model, the subsequent step involves executing IK, ID, and optimization processes to accurately compute the required muscle forces. Details of these three steps are elaborated in the subsequent subsections.

Inverse kinematics

IK analysis is used to calculate the kinematic parameters (position, velocity, and acceleration) at every time step, using the measured position of each body segment obtained from a motion capture system. While forward kinematics determines the orientation of body segments based on known joint angles, IK reverses this process to determine the joint angles required to attain a certain position [119]. IK is often achieved by reducing the disparity between the locations of the experimental and model markers for each time frame. The minimum least squares equation solved during IK is depicted in Eq. (7), where $\:i$ and $\:j$ denote the number of markers and joint angles, respectively, $\:{\varvec{x}}_{\varvec{i}}^{\varvec{e}\varvec{x}\varvec{p}}$ represents the experimental marker position, $\:{\varvec{x}}_{\varvec{i}}\left(\varvec{q}\right)$ indicates the position of the corresponding model marker $\:i$, $\:{w}_{i}$ denote marker weights, $\:\varvec{q}$ represents the generalized vector of joint angles, $\:{\varvec{q}}_{\varvec{j}}^{\varvec{e}\varvec{x}\varvec{p}}$ represents the experimental value of the joint angle, and $\:{\omega\:}_{i}$ indicates the weight of joint angles. It should be noted that $\:{\varvec{q}}_{\varvec{j}}={\varvec{q}}_{\varvec{j}}^{\varvec{e}\varvec{x}\varvec{p}}$ for all prescribed joint angles. After determining the position and joint angle of each body segment, velocity and acceleration are computed by differentiation with respect to time [120].

$$\mathop {\min }\limits_{\rm{q}} \left[ {\sum\limits_{i \in {\rm{markers}}} {\;{w_i}\parallel {x_i}^{\exp } - {x_i}(q){\parallel ^2} + \sum\limits_{j \in {\rm{joint}}\;{\rm{angles}}} {{w_j}\parallel q_j^{\exp } - {q_j}{\parallel ^2}} } } \right]$$

(7)

Inverse Dynamics

ID analysis utilizes the calculated segment kinematics (positions, velocities, accelerations), as well as the measured GRFs, to determine the muscle and joint reaction forces at each time frame. Rooted in Newton’s second law of motion, the general equation of motion governing the kinematics and kinetics of MSK models obtained through ID is depicted in Eq. (8), where $\:M$ and $\:R$ represent the muscle and joint reaction forces, respectively, while $\:\varvec{C}$ is the coefficient matrix given by $\:\varvec{C}=\left[\begin{array}{cc}{\varvec{C}}^{\left(M\right)}&\:{\varvec{C}}^{\left(R\right)}\end{array}\right]$. $\:\varvec{f}$ represents the unknown forces vector, where $\:\varvec{f}={\left[\begin{array}{cc}{\varvec{f}}^{\left(M\right)}&\:{\varvec{f}}^{\left(R\right)}\end{array}\right]}^{T}$, $\:\varvec{M}$ is the mass matrix $\:,$$\:\varvec{b}$ denotes the vector of gyroscopic terms, while $\:{\varvec{g}}^{\left(app\right)}$ represents the applied forces and moments [113, 121].

The system is statistically determinate if the force matrix has the same number of elements as in Eq. (8). However, the sophisticated intricate design of the human MSK system typically results in redundancy or statically indeterminate system, where the unknown force matrix has more elements than Eq. (8).

$$\:\left[\begin{array}{cc}{\varvec{C}}^{\left(M\right)}&\:{\varvec{C}}^{\left(R\right)}\end{array}\right]\left[\begin{array}{c}{\varvec{f}}^{\left(M\right)}\\\:{\varvec{f}}^{\left(R\right)}\end{array}\right]=\varvec{M}\dot{\varvec{v}}+\varvec{b}-{\varvec{g}}^{\left(app\right)}$$

(8)

Muscle recruitment

The human MSK system is statistically over determinate because it has more elements in the force matrix than DOF. Therefore, an infinite number of muscle and joint reaction force values can satisfy Eq. (8). This issue is typically resolved in MSK modeling by formulating an optimization problem that enables optimal muscle recruitment:

$$\min H({f^{(m)}})Cf = M\dot v + b - {g^{(app)}}s.t.0 \le {f^{(m)}} \le {s^{(M)}}$$

(9)

In this formulation, $\:H$ is a scalar objective/cost function which must be minimized while satisfying the equilibrium equation. The constraints state that the muscles only generate pulling forces, and their maximum force is constrained by their instantaneous strength ($\:{\varvec{s}}^{\left(M\right)}$), hence satisfying physiological fidelity.

Many optimization criteria for muscle recruitment have been proposed in the literature; however, the most widely used is the polynomial criterion, expressed by Eq. (10), where $\:p$ represents the power of the polynomial [113]. Description of other criteria such as soft saturation, as well as the min/max criterion, can be found in [121].

$$\:H\left({\varvec{f}}^{\left(M\right)}\right)=\sum\:_{i=1}^{{n}^{\left(M\right)}}{\left(\frac{{\varvec{f}}_{i}^{\left(M\right)}}{{\varvec{s}}_{i}^{\left(M\right)}}\right)}^{p}$$

(10)

Using the polynomial criterion as the objective function for muscle recruitment is referred to as static optimization. Despite its frequent use, this approach has been criticized on several points: 1- Kinematic Accuracy: The precise measurement and processing of segment kinematics is the first prerequisite for ensuring accuracy of ID analysis. 2- Physiological characteristics: the accuracy of muscle force prediction is determined by the underlying physiological characteristics of the muscle [122]. 3- Passive forces: this method does not consider the contribution of passive forces [123].

The CMC algorithm, developed by Thelen and Anderson [124], is another muscle recruitment technique available in OpenSim. This algorithm employs the predicted muscle excitations from static optimization in a forward simulation to predict desired accelerations for the subsequent time step [123]. The CMC algorithm, therefore, addresses some of the limitations of static optimization by incorporating dynamic elements into the muscle recruitment process. More specifically, a static optimization problem is solved to calculate the muscle excitations which yield the desired accelerations derived from IK analysis for tracking experimental motion. Importantly, the CMC algorithm integrates a proportional-derivative control law to account for the differences between experimental and model kinematics, hence ensuring more accurate motion tracking, whilst considering the delay between muscle activation and force production. This is achieved by integrating the dynamics of muscle contraction and activation, towards predicting the range of potential muscle forces that could be generated in the next time step. Unlike static optimization, which primarily considers active muscle forces, the CMC method also accounts for the contribution of passive muscle forces, leading to improved steady state muscle force predictions.

Roelker et al. [123] compared static optimization with the CMC algorithm, utilizing Gait2392 and FBRM in OpenSim. Their study found that CMC resulted in greater muscle activations, co-contraction indices, and forces as compared to static optimization. The results indicated that muscle force and activation predictions obtained by the CMC algorithm were more sensitive to model selection than those generated from static optimization.

EMG-driven and EMG-informed MSK models

EMG-driven MSK models offer a valuable alternative to optimization-based approaches for predicting muscle forces. These models utilize EMG activation data and three-dimensional joint angles to drive MSK simulations, providing subject-specific and EMG data-driven realistic predictions of muscle activity [125, 126]. The primary advantage of such models is their integration of measured EMG data, which enhances the physiological accuracy of muscle force predictions. On the other hand, they face several challenges including limited EMG data, typically measured from limited number of superficial muscles, as well as EMG data prone measurement and normalization errors due to noise and artifacts [127]. These limitations restrict the capability of EMG-driven models to predict joint moments across multiple DOF.

To address these drawbacks, EMG-informed models emerged based on integrating EMG data with traditional optimization techniques. Sartori et al. [128] proposed such model, where EMG data is used to inform and constrain the optimization process, thereby improving the accuracy of muscle force predictions. Pizzolato et al. [127] compared EMG-driven, EMG-informed, and static optimization models. They indicated that EMG-informed models provided more accurate predictions of muscle excitations and forces as compared to other methods. Similarly, Banks et al. [129] used an EMG-informed model to estimate internal lower back muscular demands accurately. While these and other pioneers have attempted to decipher and understand the complex muscle recruitment strategies used in the human body, there remains a need for a thorough evaluation of existing strategies to determine the most effective technique for muscle recruitment [113].

It should be noted that AMS uses the min/max criterion [130] for muscle recruitment, whereas both static optimization and CMC options are available in OpenSim [123]. Furthermore, OpenSim supports EMG-driven and EMG-informed techniques through the Calibrated EMG-Informed Neuromusculoskeletal Modeling (CEINMS) toolbox [127].

Post processing in gait analysis

The final phase of gait analysis using MSK modeling involves post-processing to derive the desired outcomes. This process begins with identifying the gait cycles and concludes with averaging the variables across multiple gait cycles. Figure 8 provides an overview of the potential outputs from experimental and simulated gait analysis. Notably, new additions to conventional parameters include muscle parameters and internal JRFs. The summary of studies focusing on the outputs of gait analysis using MSK modeling is presented in Supplementary Table S5. This section describes the specifics of gait cycle identification, JRFs, and muscle parameters.

Fig. 8

An overview of the parameters that could be derived from the MSK modeling-based gait analysis [119]

Gait cycle detection

Gait cycle detection using GRFs typically involves analyzing changes in the center of pressure in the medial-lateral direction. Several studies have proposed different methods for estimating GRF components during the gait cycle. A common approach is to identify specific instances of heel strike-toe-off based on the vertical component of the GRF. A heel strike is typically detected when there is a sudden increase in the vertical GRF component, usually between 10 and 20 N. Conversely, toe-off based is identified when the vertical GRF component drops below a predefined threshold [26, 131,132,133]. The precise moments within the gait cycle are critical for segmenting the data into complete cycles for further analyses.

Joint reaction forces and moments (JRF&Ms)

JRF&Ms are essential outcomes of gait analysis, providing important insight into the biomechanical forces and torques experienced by joints, such as the hip, knee, ankle, and any other joints (if specified in the MSK model) during walking. These factors are vital for comprehending muscle biomechanics and function, joint loading, as well as overall movement efficiency.

Hip Joint

The JRFs during gait at the hip joint represent the forces applied to the femoral head due to load transfer from the pelvis to the femur while walking. These forces are influenced by geometry, muscle contractions, joint kinematics, and GRFs. In addition, JRMs at the hip characterize the torques responsible for rotational movement along the joint axis. Quantifying hip JRF&Ms is critical for assessing hip joint stability, muscle coordination, and potential risk factors for hip-related injuries or diseases. For example, Skubich and Piszczatowski [134] utilized MSK modeling to calculate the hip muscle and JRFs during walking in a study involving 99 healthy participants (45 males and 54 females, ages 18 to 36) with no history of prior injuries. Their results revealed that the total, medio-lateral, and proximo-distal waveforms of the hip JRFs exhibit two peaks during the gait cycle separated by a local minimum.

Knee Joint

At the knee joint, JRF&Ms provide relevant insight into load distribution between the femur and tibia, as well as that between the femur and patella. Knee JRMs reflect the rotational torques impacting knee flexion and extension in the sagittal plane, as well as abduction/adduction moments in the frontal plane, and internal/external rotational moments in the transverse plane. Knee JRFs depict the compressive and shear forces experienced by the joint in all 3 planes. These metrics are critical for evaluating knee joint biomechanics, ligamentous stability, as well as stress and strain patterns which impact knee injuries, such as cartilage and ligament tears or osteoarthritis.

Ankle Joint

JRF&Ms at the ankle joint describe the interactions between the tibia, talus, and calcaneus bones. Notably, the complexity of the ankle anatomy has a significant influence on the function and biomechanical performance of the joint. Ankle JRFs reflect the forces transferred between the bones, influenced by the bony anatomy, GRFs, muscle and ligament forces, and foot mechanics. Meanwhile, the ankle JRMs define the torques responsible for dorsiflexion and plantarflexion movements in the sagittal plane, as well as to a smaller extent inversion/eversion in the frontal plane and abduction/adduction in the transverse plane. Understanding ankle JRF&Ms is crucial for assessing ankle joint stability, foot-ground contact, and biomechanical adaptations to different walking conditions and surfaces.

Muscle activation and forces

In gait analysis using MSK modeling, muscle parameters such as forces, activation levels, power, and muscle-tendon unit (MTU) length provide important information about the biomechanical components of human movement. The number of muscles in MSK models varies; for example, TLEM 2.0 includes 166 MTUs, while the Delp model features only 43 MTUs.

Muscle forces are typically a relevant output obtained from gait analysis. Researchers can estimate the forces generated by specific muscles during the gait cycle by combining force plate data, EMG, and MSK modeling. These force computations aid in identifying muscle synergies and the key contributing muscle groups to movement and propulsion during walking, providing insights into energy expenditure and muscle efficiency.

In addition to muscle force estimation, gait analysis provides information regarding muscle activation levels, which reflect which muscles are active and indicate the degree to which these muscles engage or contract during movement. Clinicians can analyze these activation levels to assess the timing and intensity of muscle activation, which helps them comprehend muscle coordination and synchronization during walking. Additional important muscle related metrics obtained from gait analysis include muscle power and contraction velocity. Muscle power represents the rate at which muscles perform work, whereas contraction velocity describes the rate at which muscle length changes during contraction. These measurements provide insights into muscle performance, efficiency, and fatigue resistance during walking. MTU length can also be determined using MSK simulations. Understanding how MTU length changes over the gait cycle is critical for assessing muscle flexibility, joint ROM, and biomechanical functional adaptations to various walking patterns.

Gait analysis using MSK modeling has been applied in various clinical research contexts. Some examples include investigating the relationship between forward and backward gait [135]; exploring inter-individual differences and similarities in muscle forces at the ankle joint [136]; examining the forces produced by the ankle joint’s triceps surae in forward propulsion [137]; analyzing the effects of perturbed walking on muscle forces in males and females [138]; and studying the influence of different modes of locomotion (walking and running) on the lengths of the hip MTUs [139].

MSK modeling applications in gait analysis

Over the past decade, MSK modeling has significantly advanced the understanding of both healthy and impaired human biomechanics. This powerful tool, which combines multibody dynamics calculations with imaging techniques, has proven invaluable, not only in the design and development of implants and prostheses, but also in clinical data-driven informed decision-making. Gait analysis using MSK modeling is increasingly employed to assess and understand the biomechanics for individuals with various health conditions, including post stroke movement dysfunction, Parkinson’s’ disease, cerebral palsy, diabetes, obesity, as well as ageing. The details of studies on the applications of MSK modeling in gait analysis are discussed in this section and summarized in Supplementary Tables S6 to S10.

Hip impairments

Gait analysis plays a crucial role in the research and management of hip impairments, providing insight into joint biomechanics and muscle function. For example, in the case of Acetabular dysplasia, a condition characterized by a shallow, poorly formed, or mispositioned acetabulum (hip socket) [140, 141] which significantly affects the biomechanics of the hip joint, gait analysis can quantify the JRF&Ms and muscular forces, offering valuable data on the range of motion and potential joint overloading. This information is essential for identifying areas susceptible to injury and guiding targeted treatment to prevent further joint deterioration.

Gait analysis is also instrumental in evaluating the effects of femoral offset changes in patients undergoing THA [142]. Femur offset, defined as the lateral distance from the center of the hip to the shaft of the femur, significantly influences hip stability, muscle function, and overall walking performance. By assessing abductor muscle strength, JRFs, and joint moments, gait analysis provides critical information for preoperative planning, implant selection, and postoperative rehabilitation, ensuring optimal patient outcomes. Femoro-acetabular impingement (FAI) refers to the abnormal contact between the femoral head and the acetabulum, leading to joint damage and altered gait patterns [143]. Gait analysis helps in evaluating muscle forces and joint loading patterns in individuals with FAI, revealing key changes in walking mechanics, muscle imbalances, and compensatory strategies. This information is useful for developing personalized rehabilitation programs and surgical interventions aimed at restoring normal hip function.

For patients undergoing metal-on-metal resurfacing hip arthroplasty (RHA), gait analysis aids in assessing the biomechanical performance of the hip joint [144]. By investigating the relationship between metal ion levels and hip loading patterns, researchers can predict potential implant wear and corrosion rates. This information guides implant design, material selection, implant monitoring, and patient management regimens, ensuring safety and optimal lifespan of hip implants.

Knee impairments

Gait analysis provides a detailed understanding of knee joint biomechanics by quantifying the kinematics, kinetics, and muscle activation patterns during walking. Such in-depth analysis helps elucidate the complex interplay between joint structure, function, and disease. Key applications of gait analysis in knee impairments include identification of abnormal movement patterns associated with various knee conditions and pathologies, such as osteoarthritis (OA), cartilage degeneration, ligamentous injuries, patellofemoral injuries, and post-surgical complications [65, 145,146,147,148]. By examining joint angles and moments, gait analysis can detect deviations from normal walking patterns, which are critical for diagnosing and managing these conditions. In patients with knee OA, gait analysis helps understand the impact of the disease on joint biomechanics. Such as revealing how OA affects joint loading movement efficiency. Rossom et al. [145] studied the impact of end-stage knee and hip OA on JRFs and found that both groups exhibited significantly reduced load bearing capability on the affected joint during walking. This reduction in load is often a compensatory mechanism to minimize pain and joint damage.

Similarly, Richards et al. [149] used MSK modeling to investigate the impact of varying knee OA levels of severity on predicted muscle forces and JRFs. This study demonstrated that patients with moderate knee OA exhibited lower peak knee JRF as compared to those of healthy participants. This suggests that moderate OA may not drastically alter knee mechanics as compared to severe OA, although that depends on other factors, including age. Meireles et al. [150] reported similar observations and concluded that the JRFs of the knee remain relatively unaltered during the early stages of OA. This implies that significant alterations in knee joint mechanics might only occur as OA progresses to more advanced stages. Farrokhi et al. [23] employed gait analysis to investigate the impact of intermittent versus continuous walking on knee pain and contact forces in individuals with knee OA. Their findings revealed that continuous walking increased knee pain, whereas intermittent walking did not exacerbate pain levels.

Muscle function and activation patterns are also critical components of knee impairment research, and gait analysis is a valuable tool for studying these factors. Huang et al. [151] investigated the impact of anterior cruciate ligament (ACL) deficiency and dysfunction (ACLDM) on lower body muscle activation patterns, reporting that patients with ACL rupture and ACLDM exhibited reduced gastrocnemius muscle force as compared to healthy subjects. On the other hand, Aghdam et al. [152] reported that patients with an ACL rupture experienced greater knee joint loading than healthy individuals. Furthermore, Kang et al. [153] investigated the effects of posterior cruciate ligament (PCL) deficiency on knee JRFs, and found that knee JRFs increased during gait as the PCL’s deficiency worsened.

Therapeutic interventions and rehabilitation represent another significant area where gait analysis is extensively utilized in knee associated impairment. Schroeder et al. [68] employed gait analysis to investigate muscle loads during gait post-ACL surgery, revealing that knee muscle forces significantly decreased following ACL reconstruction. Similarly, Richards et al. [154] examined the relationship between knee JRF&Ms in individuals with knee OA during normal and modified gait (wide-steps and toe-in). Their study concluded that while both wide-steps and toe-in gait modifications did not decrease the overall JRFs of the knee as compared to normal walking, they significantly reduced the ratio of medial and total knee JRFs. Additionally, Thorsen et al. [155] explored the variations in knee muscle forces and mediolateral tibiofemoral compressive forces among patients with TKA during gait. They found no statistical differences in muscle or compressive forces between limbs.

In summary, gait analysis offers a comprehensive, quantitative, and multidimensional approach to studying knee-associated impairments by providing critical biomechanical data essential for gait and motion dysfunction assessment treatment planning, assessment of therapeutic efficacy, as well as advancement in biomechanical research.

Foot impairments

Gait analysis plays a crucial role in foot impairment research by evaluating ankle muscle forces and their impact on stress distribution within the foot and ankle joints. Hejazi and Rasmussen [156] used finite element analysis and MSK modeling to investigate the impact of ankle muscle forces on stress distribution in patients with ankle arthrodesis. Their findings indicated that arthrodesis patients experienced reduced stress on the medial malleolus when the gastrocnemius-soleus muscle activity increased, as compared to healthy individuals. Moreover, gait analysis enables the comparison of pelvic kinematics and lumbopelvic muscle forces between individuals with foot hyperpronation and those with properly aligned feet. Yazdani et al. [157, 158] investigated these differences, revealing the altered biomechanical and muscle activation patterns in individuals with foot hyperpronation. Their research revealed that subjects with hyperpronation exhibited a greater anterior pelvic tilt throughout 20–80% of the gait stance phase as compared to subjects with correctly aligned feet.

In individuals with chronic ankle instability, gait analysis can be used to investigate the effects of vibration-based treatment on ankle JRFs and vertical GRFs during walking. Studies by Jang et al. [159] and Jang and Wikstrom [160] demonstrated that ankle JRFs and vertical GRFs dramatically reduced following vibration-based treatment. On the other hand, individuals with chronic ankle instability exhibited statistically significant larger anteroposterior and mediolateral shear forces (initial peak and impulse) as compared to healthy control, while their compressive forces were lower.

Greve et al. [161] evaluated the impact of rocker shoes on plantar aponeurosis strain in patients with plantar fasciitis as compared to a healthy control group during locomotion. The study found that wearing rocker shoes did not result in a statistically significant reduction in plantar aponeurosis strain.

Gait analysis can also be used to investigate how weak ankle dorsiflexion muscles affect the forces exerted by other muscle groups during walking. Błażkiewicz and Wit [14] conducted a study to examine these relationships and observed that the increased muscle forces of the flexor hallucis, tibialis posterior, and flexor digitorum longus on the ankle joint, as well as the semimembranosus and semitendinosus at the surrounding joints, compensated for the weakness of the ankle dorsiflexion muscle.

Orthosis and surgical procedures evaluation

MSK modeling for gait analysis extends its utility to evaluating the effectiveness of various medical protocols and interventions, including arthroplasty surgery techniques [162], functional resistance training [163], knee orthotics and braces [164], footwear [165], spinal devices and procedures [166], and foot orthoses [167].

Many researchers employed MSK modeling during gait towards the evaluation of footwear interventions and their impact on GRFs, spatiotemporal features, and stability in specific populations. Li et al. [165] investigated the impact of normal, negative, and positive heel shoes on stride biomechanics in pregnant women. They found that while the center of pressure reduced dramatically with negative heel shoes, step length and walking velocity increased with heel-toe drop. Additionally, negative and positive heel shoes exhibited significantly higher and lower vertical GRFs, respectively, as compared to regular shoes.

Wang et al. [166] studied the impact of spinal orthosis on trunk and lower extremity joint loads using MSK modeling during gait. They found that the maximum compressive forces and moments at the hip, knee, and ankle joints did not significantly differ between conditions of wearing spinal orthosis, assuming normal posture, or adopting a poor posture. Likewise, MSK modeling offers a valuable tool for investigating the effects of knee braces and the perturbations they generate on muscle forces and joint mechanics during walking. Yap et al. [164] investigated the impact of perturbations generated by knee braces on knee muscle forces. They found that normal knee exhibited higher forces in the knee flexor muscles than the braced knee. In contrast, the braced knee demonstrated higher forces in the knee extensor muscles as compared to the normal knee.

Marconi et al. [167] employed MSK simulations to investigate the effect of mass distribution of powered ankle foot orthosis on lower body muscle forces and joint moments. Their findings revealed that while the additional mass of the orthoses did not notably alter ankle joint moment, it significantly altered hip and knee joint moments throughout 60–100% of the gait cycle. Additionally, MSK modeling proves to be a valuable tool for evaluating the effects of functional resistance training, as demonstrated by Washabaugh et al. [163]. Their study investigated how functional resistance training altered lower-limb muscles and JRFs during walking. Notably, they observed that ankle weights lead to a substantial increase in knee flexion and hip extension moments by 174% and 477%, respectively.

MSK modeling has also been instrumental in exploring the biomechanical alterations and functional adaptations associated with various hip arthroplasty surgical methods. Wesseling et al. [162] investigated three patient groups representing distinct hip arthroplasty surgery techniques and revealed that, irrespective of the surgical approach, individuals with hip impairments exhibited reduced hip JRFs during walking as compared to their healthy counterparts.

Other applications of MSK modeling for gait analysis

Beyond analyzing lower limb (hip, knee, and ankle) biomechanical function, gait analysis has emerged as a diverse and valuable method for investigating the impact of various medical conditions on gait. These conditions include cerebral palsy [168,169,170,171], diabetes [172, 173], obesity [174], aging [175], lumbar disc herniation [176], stroke [177], Osgood-Schlatter disease [178], and Parkinson’s disease [21, 22, 179].

Salami et al. [168] compared trunk kinematics and hip abduction moments in children with cerebral palsy to healthy controls. They found that the lumbar ROM in the cerebral palsy group was double that of the healthy control group.

In studying the impact of obesity on joint forces during walking, Haight et al. [174] investigated tibiofemoral joint forces during fast-level and slow-uphill walking. They observed that obese subjects exhibited higher compressive tibiofemoral forces as compared to non-obese subjects in both walking conditions. Khalaf et al. [180] investigated the influence of overweight (BMI ≥ 25) on the dynamic peak plantar pressure (PPP) distribution during gait using traditional and nonlinear dynamic measures in young college students. This study found that the overweight group had significantly higher mean PPP and entropy values as compared to the control group. The overweight group also had significantly higher PPP underfoot as compared to heel. The authors concluded that the overweight group was at risk of ulceration and that different shoe wear should be considered. Gait analysis also facilitates the investigation of age-related differences in muscle forces during walking. Toda et al. [175] reported that muscle coordination around the ankle changes with age, with lower limb muscles co-contraction occurring during the late stance phase.

Gomes et al. [172] used MSK modeling to calculate muscular forces in diabetic individuals with and without polyneuropathy. Their findings indicated that diabetic patients without polyneuropathy exhibited decreased peak soleus muscle force during middle and late stance phases as compared to healthy individuals. Conversely, diabetic patients with polyneuropathy showed a decreased gastrocnemius medialis muscle forces as compared to those without polyneuropathy. Similarly, Scarton et al. [173] examined the impact of diabetic peripheral neuropathy on muscle force estimation in the lower body using gait analysis. This study revealed statically significant differences in muscle forces, joint kinematics, and JRFs between patients with diabetic peripheral neuropathy and a healthy control group. Notably, the patient group demonstrated weaker ankle muscle forces than the healthy controls. Khalaf et al. [181] studied gait alterations in the UAE population with and without diabetic complications using both traditional and entropy measures. Their results indicated that intrinsic sensorimotor feedback plays a vital role in regulating gait for patients with peripheral neuropathy.

Gait analysis using MSK modeling also enables the estimation of lower body power and muscular forces in individuals with peripheral arterial disease. Rahman et al. [182] investigated the effect of peripheral artery disease on gait and discovered that the impaired group’s knee extensors, plantar flexors, and hip flexors exhibited significantly lower power and force as compared to healthy individuals. Furthermore, gait analysis using MSK modeling serves as a valuable tool for assessing the impact of minimally invasive surgery on spatiotemporal characteristics and gait kinetics in individuals with lumbar disc herniation. Wang et al. [176] used MSK modeling to investigate these relationships, revealing lower peak muscle forces of the rectus femoris, tensor fasciae latae, and biceps femoris long-head on the affected side, as well as reductions in gait cycle time and stride length post-surgery.

In summary, gait analysis using MSK modeling offers a comprehensive, quantitative, and multidimensional approach to examining various biomechanical challenges and health conditions, providing data-driven relevant insights for disease characterization, as well as rehabilitation and treatment optimization.

Discussion

The integration of multibody dynamics with biofidelic MSK modeling allows for significant advancement in motion and gait assessment, offering a valuable alternative to traditional methods, which often rely on complex costly instrumentation and/or invasive procedures. A notable advantage to this approach lies in its ability to predict internal JRF&Ms and muscle forces using kinematics alone. This enables the quantification of internal parameters which cannot be measured non-invasively in-vivo thereby expanding the scope and accessibility of gait analysis in both research and clinical settings.

Importantly, MSK modeling in gait analysis has proven to be a versatile and effective method for studying a wide range of functional biomechanical challenges and health conditions. The investigations explored in this systematic review reveal the successful quantification of parameters non-invasively and beyond what is possible with in-vivo measurements. The development, enhancement, and clinical applicability of such in-silico MSK based tools hold significant potential to enhance patient care and open doors to a paradigm shift in healthcare towards personalized, precision, preventative and predictive medicine. This can improve the entire healthcare cycle from data-driven diagnostics to informed decision-making and quantitative evaluation of interventions.

This review aimed to provide readers with a concise overview of the latest state-of-the-art research in the development and application of MSK models. In biomechanical research, MSK modeling has emerged as an invaluable tool for addressing “what if” type of questions in a cost and time effective manner, providing researchers and clinicians with means for investigating various treatment options and scenarios towards optimal ones. On the other hand, it is crucial to critically evaluate the limitations inherent to current MSK models to ensure adequate accuracy, reliability, and applicability of research findings. This section, therefore, highlights the key limitations of rigid body musculoskeletal models and outlines potential areas for future research.

Limitations of MSK modeling in clinical practice

Despite the widespread use of MSK modeling in research, the transition of MSK modeling into clinical practice has been limited by several critical challenges. The first and most significant limitation is the dependence on bulky and cumbersome equipment for collecting kinematic data, which serves as the input for inverse dynamic simulations. For instance, reflective markers or IMU sensors must be placed on subjects, a process that is often uncomfortable and challenging to execute. Additionally, the use of force plates to capture GRFs presents another hurdle. Force plates can inadvertently alter a subject’s natural gait or limit data collection to just a few steps. While artificial intelligence offers potential solutions for predicting kinematics and GRFs, this technology is still immature and requires further development and research to become a viable alternative.

One of the primary obstacles is the difficulty in achieving comprehensive personalization of MSK models to reflect the unique anatomical, physiological, and neurological characteristics of individual patients. The inherent variability in human anatomy and physiology, influenced by factors such as the specific cause and type of impairment, can significantly affect the outcomes of rehabilitation, treatment, or surgical interventions. Traditional scaling methods, often based on generic cadaveric models, may not sufficiently capture these individual differences. Although advanced techniques, such as those utilizing MRI or CT-based scaling, provide a higher degree of personalization, they are often labour-intensive, time-consuming, and costly. This limitation underscores the necessity for the development of more efficient, cost-effective, and user-friendly scaling methodologies that can streamline the personalization of MSK models, thereby enhancing their applicability in clinical settings.

Another significant limitation is the challenge of predicting post-treatment function, given the complex interplay between human movement and the neuromuscular system. Human movement is not purely a mechanical process; it is dynamically modulated by neuromuscular control, which can adapt over time in response to treatment or disease progression. However, most existing MSK models do not account for this functional adaptability, which restricts their ability to predict long-term outcomes or accommodate changes in patient condition. Incorporating neuromuscular adaptability into MSK models is therefore critical for making these models more reflective of real-world scenarios [183].

A salient limitation of current MSK modeling frameworks arises from the foundational concepts which underlie many models available in prominent software platforms, such as OpenSim and AMS. OpenSim models prioritize inclusivity and representativeness, drawing upon data from a diverse range of subjects to encapsulate broader population dynamics. By contrast, AMS models emphasize data consistency, advocating for the integration of cohesive datasets rather than combining fragmented data from disparate studies [74, 75]. This divergence in conceptual frameworks introduces inherent biases and limitations that researchers must carefully consider when selecting and deploying MSK models for biomechanical simulations.

Furthermore, the current muscle recruitment criteria employed in MSK models warrant critical reevaluation and methodological refinement. While static optimization remains a prevalent approach, its dependency on precise kinematic measurements and its exclusion of passive muscle forces compromises its comprehensiveness and predictive accuracy. Alternative strategies, such as EMG-driven or EMG-informed techniques, offer viable solutions but come with their own set of challenges. EMG-driven models, for example, often impose substantial computational burdens and may not be applicable across all muscle groups due to limitations in measurable EMG data and potential normalization errors. Optimizing muscle recruitment strategies is therefore in need of rigorous exploration. This entails conducting thorough methodological comparisons and validations to determine the most suitable and reliable technique for accurately simulating muscle activations and forces within MSK models.

Additionally, the validity of MSK models remains a pressing concern. The accuracy of predictions regarding joint reaction forces and muscle forces is often limited by the scarcity of in vivo data for validation. While instrumented implants provide some insights into internal forces, the data are typically confined to a small cohort of impaired individuals and particular age groups, which limits their generalizability. The lack of in vivo muscle force data further compromises the predictive accuracy of MSK models, which is essential for reliable clinical decision-making.

In summary, while MSK modeling holds significant potential for enhancing patient care, its clinical application remains constrained by challenges related to model personalization, adaptability, and validation. Addressing these limitations is crucial for advancing the use of MSK models in personalized and effective patient treatment strategies. Outlook and Future Work

The use of MSK models in gait analysis has greatly advanced the understanding of human biomechanics and patterns of motion. However, a notable limitation of current MSK models, particularly those available in software platforms like OpenSim and AMS, is their reliance on male anatomy or geometry. Despite this, about 57.75% of the reviewed studies (67 out of 116) used these male-based models to simulate female gait and locomotion. This raises a critical question regarding the necessity and relevance of the development of a female specific MSK model.

To address the necessity of developing female specific MSK models, it is crucial to consider the physical differences between males and females and their potential impact on joint kinetics during walking. Fischer and Mitteroecker [184] studied sexual dimorphism and allometry in the human pelvis, examining a sample of 99 participants, including 46 males and 53 females. This geometric morphometric study indicated considerable anatomical differences between the male and female pelvis. The study revealed significant anatomical variations between genders. Specifically, females had shorter and more laterally extending iliac blades, a rounder and outward-projecting pelvic inlet, a shorter sacrum, and a wider subpubic angle than males. These morphological differences are substantial, with the sexual dimorphism vectors for both genders showing significant discrepancies equivalent to four standard deviations. These findings underscore substantial morphological discrepancies between the male and female pelvis, which likely affect gait mechanics and joint kinetics.

Delving deeper into the biomechanical implications of these anatomical variations, Lewis et al. [157] explored their impact on the joint ROM during gait. Their study, encompassing 44 participants (22 males and 22 females), revealed notable differences in pelvic positioning and movement between genders. On average, females exhibited a more anteriorly tilted pelvis, approximately 4 degrees more in the sagittal plane, as compared to the near-neutral pelvic position observed in males. Additionally, females demonstrated greater pelvic excursion in the frontal plane, with a statistically significant mean difference of 1.9 degrees. In the transverse plane, females also displayed enhanced pelvic mobility compared to their male counterparts. Building upon these insights, Ismail and Lewis [185] used MSK simulations to investigate the influence of pelvic tilt variances between genders on hip JRFs during gait. Their findings revealed that pelvic tilt exerted a substantial influence on hip JRFs.

In light of these findings, the absence of a female specific MSK model emerges as a notable limitation in accurately predicting JRFs and understanding biomechanical adaptations in female populations. Addressing this gap is crucial for enhancing the precision and applicability of MSK simulations in diverse research and clinical contexts.

Conclusion

Multibody MSK modeling stands as a significant advancement in the quantification of human motion, especially in gait analysis for data-driven diagnostics, medical device development and informed clinical decision-making. However, to fully realize its potential, further research is necessary to enhance the accuracy and biofidelity of these models. This includes the development of female population-specific MSK models to better represent the biomechanics of female individuals. The creation of subject-specific scaling tools which account for individual anthropometrics is also important towards more accurate personalized and precise simulations.