Inter-individual variation in vertebral kinematics affects predictions of neck musculoskeletal models

Experimental studies have found significant variation in cervical intervertebral kinematics (IVK) among healthy subjects, but the effect of this variation on biomechanical properties, such as neck strength, has not been explored. The goal of this study was to quantify variation in model predictions of extension strength, flexion strength and gravitational demand (the ratio of gravitational load from the weight of the head to neck muscle extension strength), due to inter-subject variation in IVK. IVK were measured from sagittal radiographs of 24 subjects (14F, 10M) in five postures: maximal extension, mid-extension, neutral, mid-flexion, and maximal flexion. IVK were defined by the position (anterior-posterior and superior-inferior) of each cervical vertebra with respect to T1 and its angle with respect to horizontal, and fit with a cubic polynomial over the range of motion. The IVK of each subject were scaled and incorporated into musculoskeletal models to create models that were identical in muscle force- and moment-generating properties but had subject-specific kinematics. The effect of inter-subject variation in IVK was quantified using the coefficient of variation (COV), the ratio of the standard deviation to the mean. COV of extension strength ranged from 8% to 15% over the range of motion, but COV of flexion strength was 20–80%. Moreover, the COV of gravitational demand was 80–90%, because the gravitational demand is affected by head position as well as neck strength. These results indicate that including inter-individual variation in models is important for evaluating neck musculoskeletal biomechanical properties.     

Approximating complex musculoskeletal biomechanics using multidimensional autogenerating polynomials

Computational models of the musculoskeletal system are scientific tools used to study human movement, quantify the effects of injury and disease, plan surgical interventions, or control realistic high-dimensional articulated prosthetic limbs. If the models are sufficiently accurate, they may embed complex relationships within the sensorimotor system. These potential benefits are limited by the challenge of implementing fast and accurate musculoskeletal computations. A typical hand muscle spans over 3 degrees of freedom (DOF), wrapping over complex geometrical constraints that change its moment arms and lead to complex posture-dependent variation in torque generation. Here, we report a method to accurately and efficiently calculate musculotendon length and moment arms across all physiological postures of the forearm muscles that actuate the hand and wrist. Then, we use this model to test the hypothesis that the functional similarities of muscle actions are embedded in muscle structure. The posture dependent muscle geometry, moment arms and lengths of modeled muscles were captured using autogenerating polynomials that expanded their optimal selection of terms using information measurements. The iterative process approximated 33 musculotendon actuators, each spanning up to 6 DOFs in an 18 DOF model of the human arm and hand, defined over the full physiological range of motion. Using these polynomials, the entire forearm anatomy could be computed in <10 μs, which is far better than what is required for real-time performance, and with low errors in moment arms (below 5%) and lengths (below 0.4%). Moreover, we demonstrate that the number of elements in these autogenerating polynomials does not increase exponentially with increasing muscle complexity; complexity increases linearly instead. Dimensionality reduction using the polynomial terms alone resulted in clusters comprised of muscles with similar functions, indicating the high accuracy of approximating models. We propose that this novel method of describing musculoskeletal biomechanics might further improve the applications of detailed and scalable models to describe human movement.     

Predicting muscle forces in gait from EMG signals and musculotendon kinematics

An EMG-driven muscle model for determining muscle force-time histories during gait is presented. The model, based on Hill's equation (1938), incorporates morphological data and accounts for changes in musculotendon length, velocity, and the level of muscle excitation for both concentric and eccentric contractions. Musculotendon kinematics were calculated using three-dimensional cinematography with a model of the musculoskeletal system. Muscle force-length-EMG relations were established from slow isokinetic calibrations. Walking muscle force-time histories were determined for two subjects. Joint moments calculated from the predicted muscle forces were compared with moments calculated using a linked segment, inverse dynamics approach. Moment curve correlations ranged from r = 0.72 to R = 0.97 and the root mean square (RMS) differences were from 10 to 20 Nm. Expressed as a relative RMS, the moment differences ranged from a low of 23% at the ankle to a high of 72% at the hip. No single reason for the differences between the two moment curves could be identified. Possible explanations discussed include the linear EMG-to-force assumption and how well the EMG-to-force calibration represented excitation for the whole muscle during gait, assumptions incorporated in the muscle modeling procedure, and errors inherent in validating joint moments predicted from the model to moments calculated using linked segment, inverse dynamics. The closeness with which the joint moment curves matched in the present study supports using the modeling approach proposed to determine muscle forces in gait.     

Estimation of muscle response using three-dimensional musculoskeletal models before impact situation: a simulation study

When car crash experiments are performed using cadavers or dummies, the active muscles' reaction on crash situations cannot be observed. The aim of this study is to estimate muscles' response of the major muscle groups using three-dimensional musculoskeletal model by dynamic simulations of low-speed sled-impact. The three-dimensional musculoskeletal models of eight subjects were developed, including 241 degrees of freedom and 86 muscles. The muscle parameters considering limb lengths and the force-generating properties of the muscles were redefined by optimization to fit for each subject. Kinematic data and external forces measured by motion tracking system and dynamometer were then input as boundary conditions. Through a least-squares optimization algorithm, active muscles' responses were calculated during inverse dynamic analysis tracking the motion of each subject. Electromyography for major muscles at elbow, knee, and ankle joints was measured to validate each model. For low-speed sled-impact crash, experiment and simulation with optimized and unoptimized muscle parameters were performed at 9.4 m/h and 10 m/h and muscle activities were compared among them. The muscle activities with optimized parameters were closer to experimental measurements than the results without optimization. In addition, the extensor muscle activities at knee, ankle, and elbow joint were found considerably at impact time, unlike previous studies using cadaver or dummies. This study demonstrated the need to optimize the muscle parameters to predict impact situation correctly in computational studies using musculoskeletal models. And to improve accuracy of analysis for car crash injury using humanlike dummies, muscle reflex function, major extensor muscles' response at elbow, knee, and ankle joints, should be considered.     

Optimising muscle parameters in musculoskeletal models using Monte Carlo simulation

The use of musculoskeletal simulation software has become a useful tool for modelling joint and muscle forces during human activity, including in reduced gravity because direct experimentation is difficult. Knowledge of muscle and joint loads can better inform the design of exercise protocols and exercise countermeasure equipment. In this study, the LifeModeler? (San Clemente, CA, USA) biomechanics simulation software was used to model a squat exercise. The initial model using default parameters yielded physiologically reasonable hip-joint forces but no activation was predicted in some large muscles such as rectus femoris, which have been shown to be active in 1-g performance of the activity. Parametric testing was conducted using Monte Carlo methods and combinatorial reduction to find a muscle parameter set that more closely matched physiologically observed activation patterns during the squat exercise. The rectus femoris was predicted to peak at 60.1% activation in the same test case compared to 19.2% activation using default parameters. These results indicate the critical role that muscle parameters play in joint force estimation and the need for exploration of the solution space to achieve physiologically realistic muscle activation.     

Muscle mass in musculoskeletal models

Most current models of musculoskeletal dynamics lump a muscle's mass with its body segment, and then simulate the dynamics of these body segments connected by joints. As shown here, this popular approach leads to errors in the system's inertia matrix and hence in all aspects of the dynamics. Two simplified mathematical models were created to capture the relevant features of monoarticular and biarticular muscles, and the errors were analyzed. The models were also applied to two physiological examples: the triceps surae muscles that plantar flex the human ankle and the biceps femoris posterior muscle of the rat hind limb. The analysis of errors due to lumping showed that these errors can be large. Although the errors can be reduced in some postures, they cannot be easily eliminated in models that use segment lumping. Some options for addressing these errors are discussed.     

Minimal models of multidimensional computations

The multidimensional computations performed by many biological systems are often characterized with limited information about the correlations between inputs and outputs. Given this limitation, our approach is to construct the maximum noise entropy response function of the system, leading to a closed-form and minimally biased model consistent with a given set of constraints on the input/output moments; the result is equivalent to conditional random field models from machine learning. For systems with binary outputs, such as neurons encoding sensory stimuli, the maximum noise entropy models are logistic functions whose arguments depend on the constraints. A constraint on the average output turns the binary maximum noise entropy models into minimum mutual information models, allowing for the calculation of the information content of the constraints and an information theoretic characterization of the system's computations. We use this approach to analyze the nonlinear input/output functions in macaque retina and thalamus; although these systems have been previously shown to be responsive to two input dimensions, the functional form of the response function in this reduced space had not been unambiguously identified. A second order model based on the logistic function is found to be both necessary and sufficient to accurately describe the neural responses to naturalistic stimuli, accounting for an average of 93% of the mutual information with a small number of parameters. Thus, despite the fact that the stimulus is highly non-Gaussian, the vast majority of the information in the neural responses is related to first and second order correlations. Our results suggest a principled and unbiased way to model multidimensional computations and determine the statistics of the inputs that are being encoded in the outputs.     

The effects of anatomical errors on shoulder kinematics computed using multi-body models

Joint motion calculated using multi-body models and inverse kinematics presents many advantages over direct marker-based calculations. However, the sensitivity of the computed kinematics is known to be partly caused by the model and could also be influenced by the participants’ anthropometry and sex. This study aimed to compare kinematics computed from an anatomical shoulder model based on medical images against a scaled-generic model and quantify the effects of anatomical errors and participants’ anthropometry on the calculated joint angles. Twelve participants have had planar shoulder movements experimentally captured in a motion lab, and their shoulder anatomy imaged using an MRI scanner. A shoulder multi-body dynamics model was developed for each participant, using both an image-based approach and a scaled-generic approach. Inverse kinematics have been performed using the two different modelling procedures and the three different experimental motions. Results have been compared using Bland–Altman analysis of agreement and further analysed using multi-linear regressions. Kinematics computed via an anatomical and a scaled-generic shoulder models differed in average from 3.2 to 5.4 degrees depending on the task. The MRI-based model presented smaller limits of agreement to direct kinematics than the scaled-generic model. Finally, the regression model predictors, including anatomical errors, sex, and BMI of the participant, explained from 41 to 80% of the kinematic variability between model types with respect to the task. This study highlighted the consequences of modelling precision, quantified the effects of anatomical errors on the shoulder kinematics, and showed that participants' anthropometry and sex could indirectly affect kinematic outcomes.

     

Validation of subject-specific musculoskeletal models using the anatomical reachable 3-D workspace

A novel metric for the validation of musculoskeletal models is proposed, the reachable 3-D workspace (RWS). This new metric was used to compare a generic model scaled in a standard manner to a more subject-specific model. An experimental protocol for assessing the RWS was performed by ten participants for four distinct hand-payload cases. In addition, isometric individual strength measurements were collected for 12 different directions. The strength of subject-specific musculoskeletal models was then computed using the following assumptions: (1) standard routines including the length-mass-fat (LMF) scaling law; (2) the isometric strengths of the muscle elements were optimized to the individual strength measurements using joint strength factors (JSF). The RWS of each participant was subsequently estimated from each of the scaling approaches, LMF and JSF, for the four load cases. The experimental RWS showed that the volume and shape decreased with increasing hand-payload for every participant. The lateral and frontal far-from-torso aspects of the RWS were reduced the most. These trends were reproduced by both strength scaling approaches, but the LMF-scaled models were not able to track the overall RWS volume decrease with increasing payload, since they proved to be weaker than the participants. On the other hand, the optimised JSF subject-specific models performed better on the prediction of the RWS for all payload cases across participants. The RWS can potentially be further used as a subject-specific musculoskeletal model validation, enabling quantification of the volume and shape differences between experimentally and model-predicted RWSs.     

Estimation of parameters in large offspring number models and ratios of coalescence times

The ratio of singletons to the total number of segregating sites is used to estimate a reproduction parameter in a population model of large offspring numbers without having to jointly estimate the mutation rate. For neutral genetic variation, the ratio of singletons to the total number of segregating sites is equivalent to the ratio of total length of external branches to the total length of the gene genealogy. A multinomial maximum likelihood method that takes into account more frequency classes than just the singletons is developed to estimate the parameter of another large offspring number model. The performance of these methods with regard to sample size, mutation rate, and bias, is investigated by simulation. The expected value of the ratio of the total length of external branches to the total length of the whole tree is, using simulation, shown to decrease for the Kingman coalescent as sample size increases, but can increase or decrease, depending on parameter values, for Λ coalescents. Considering ratios of tree statistics, as opposed to considering lengths of various subtrees separately, can yield better insight into the dynamics of gene genealogies.     

Estimation of multivariate frailty models using penalized partial likelihood

There exists a growing literature on the estimation of gamma distributed multiplicative shared frailty models. There is, however, often a need to model more complicated frailty structures, but attempts to extend gamma frailties run into complications. Motivated by hip replacement data with a more complicated dependence structure, we propose a model based on multiplicative frailties with a multivariate log-normal joint distribution. We give a justification and an estimation procedure for this generally structured frailty model, which is a generalization of the one presented by McGilchrist (1993, Biometrics 49, 221-225). The estimation is based on Laplace approximation of the likelihood function. This leads to estimating equations based on a penalized fixed effects partial likelihood, where the marginal distribution of the frailty terms determines the penalty term. The tuning parameters of the penalty function, i.e., the frailty variances, are estimated by maximizing an approximate profile likelihood. The performance of the approximation is evaluated by simulation, and the frailty model is fitted to the hip replacement data.     

Solution of large compartmental models using numerical transform inversion

Compartmental models of biological or physical systems are often described by a system of “stiff” differential equations. In this paper an algorithm for solving a system with linear coefficients is presented that employs numerical inversion of the Laplace transform of the model equations. The inversion algorithms and Gear's backward differentiation method are compared for two stiff test problems and a differential system governing a 27-compartment model of bile acid transport and metabolism. The inversion algorithm is reliable, requires modest computation time on a desktop computer and provides better accuracy than Gear's method, especially for the extremely stiff example.     

Computational method for muscle-path representation in musculoskeletal models

 This paper presents a new and efficient method to calculate the line-of-action of a muscle as it wraps over bones and other tissues on its way from origin to insertion. The muscle is assumed to be a one-dimensional, massless, taut string, and the surfaces of bones that the muscle may wrap around are approximated by cross-sectional boundaries obtained by slicing geometrical models of bones. Each cross-sectional boundary is approximated by a series of connected line segments. Thus, the muscle path to be calculated is piecewise linear with vertices being the contact points on the cross-sectional boundaries of the bones. Any level of geometric accuracy can be obtained by increasing the number of cross sections and the number of line segments in each cross section. The algorithm is computationally efficient even for large numbers of cross sections. Received: 9 July 2001 / Accepted in revised form: 11 March 2002     

Extended foot-ankle musculoskeletal models for application in movement analysis

Multibody simulations of human motion require representative models of the anatomical structures. A model that captures the complexity of the foot is still lacking. In the present work, two detailed 3D multibody foot-ankle models generated based on CT scans using a semi-automatic tool are described. The proposed models consists of five rigid segments (talus, calcaneus, midfoot, forefoot and toes), connected by five joints (ankle, subtalar, midtarsal, tarsometatarsal and metatarsophalangeal), one with 15DOF and the other with 8DOF. The calculated kinematics of both models were evaluated using gait trials and compared against literature, both presenting realistic results. An inverse dynamic analysis was performed for the 8DOF model, again presenting feasible dynamic results.     

Estimation of temporal gait parameters using Bayesian models on acceleration signals

The purpose of this study is to develop a system capable of performing calculation of temporal gait parameters using two low-cost wireless accelerometers and artificial intelligence-based techniques as part of a larger research project for conducting human gait analysis. Ten healthy subjects of different ages participated in this study and performed controlled walking tests. Two wireless accelerometers were placed on their ankles. Raw acceleration signals were processed in order to obtain gait patterns from characteristic peaks related to steps. A Bayesian model was implemented to classify the characteristic peaks into steps or nonsteps. The acceleration signals were segmented based on gait events, such as heel strike and toe-off, of actual steps. Temporal gait parameters, such as cadence, ambulation time, step time, gait cycle time, stance and swing phase time, simple and double support time, were estimated from segmented acceleration signals. Gait data-sets were divided into two groups of ages to test Bayesian models in order to classify the characteristic peaks. The mean error obtained from calculating the temporal gait parameters was 4.6%. Bayesian models are useful techniques that can be applied to classification of gait data of subjects at different ages with promising results     

A note on multidimensional models of neutral mutation

A generalisation of the Ohta-Kimura charge-state model for neutral mutation, in which alleles are represented by points in d-dimensional space, seems appropriate in the light of recent experimental developments. It is noted that existing methods of analysis extend almost trivially to this more general model, and the point is illustrated by calculating the effective number of alleles.