INVESTIGATION OF THE LOADING AT THE KNEE JOINT COMPLEX USING AN EMG-BASED CONSTITUTIVE LAW FOR SKELETAL MUSCLE FORCE
Abstract
This study investigated the predictive ability of the skeletal muscle force model presented by Knodel et al. [Knodel NB, Lawson LB, Nauman EA, “An emg-based constitutive law for force generation in skeletal muscle-part i: Model development,” J Biomech Eng (in press), doi: 10.1115/1.4053568] on the knee joint. It has previously been validated on the ankle joint [Knodel NB, Calvert LB, Bywater EA, Lamia JP, Patel SN, Nauman EA, “An emg-based constitutive law for force generation in skeletal muscle-part ii: Model validation on the ankle joint complex,” Submitted for Publication] and this paper aimed to identify how well it, and the solution process, performed on a more complex articulation. The knee joint’s surrounding musculoskeletal tissue loading was also identified. Ten subjects (five male and five female) performed six exercises targeting the muscles that cross the knee joint. Motion capture, electromyography, and force plate data was collected during the exercises for use in the analysis program written in MATLAB and magnetic resonance images were used to observe subject-specific ligament and tendon data at the knee articulation. OpenSim [Delp, SL, Anderson FC, Arnold AS, Loan P, Habib A, John CT, Guendelman E, Thelen DG, “Opensim: Open-source software to create and analyze dynamic simulations of movement,” IEEE Trans Biomed Eng 54(11):1940–1950, 2007, doi: 10.1109/TBME.2007.901024] was used for scaling a generic lower extremity anatomical model of each subject. Five of the six exercises were used to calculate each muscle’s constant, Km [Knodel NB, Lawson LB, Nauman EA, “An emg-based constitutive law for force generation in skeletal muscle-part i: Model development,” J Biomech Eng (in press), doi: 10.1115/1.4053568; Knodel NB, Calvert LB, Bywater EA, Lamia JP, Patel SN, Nauman EA, “An emg-based constitutive law for force generation in skeletal muscle-part ii: Model validation on the ankle joint complex,” Submitted for Publication], and the sixth was used as a testing set to identify the model’s predictive ability. Average percent errors ranged from 9.4% to 26.5% and the average across all subjects was 20.6%. The solution process produced physiologically relevant muscle forces and the surrounding tissue loading behaved as expected between the various exercises without approaching respective tensile strength values.
1. Introduction
Knee injuries are a common occurrence in the United States and place a substantial burden on the healthcare system.4,5,6,7 There are estimated to be 200,00–250,000 anterior cruciate ligament (ACL) injuries alone each year reaching $2 billion in medical costs.6,7 Understanding the loading conditions in the musculoskeletal tissues of the knee is critical for advancements in injury prevention and rehabilitation, especially in athletics where the knee is subjected to rigorous loading conditions.8,9,10,11 In fact, sports-related injuries accounted for half of the entire 6.5 million reported injuries summarized by Gage et al.7 To identify these loads, muscle forces at the knee must first be determined as they are the primary loading mechanism on the musculoskeletal system. The model by Knodel et al.1,2 has shown promising results at the ankle joint but has yet to be evaluated at the more complex knee joint where a larger number of muscles contribute to the joint loading and the distribution of the joint contact loading (JCL) is more difficult to model due to the complex anatomy. Identifying the predictive ability of the model at this more complex joint by observing activities that require activation of the majority of the primary muscles crossing the knee as well as co-contraction, such as squats, is an important step in detailing the limits of its application and accuracy.
The muscle force model presented by Knodel et al.1 demonstrated the ability to predict the components of →B, a vector comprised of all the directly measured terms from the second Newton–Euler equation (2) by subtracting the joint torques due to the ground reaction forces and system mass from the kinematic terms, with less than 10% error.2 This provided confidence in the skeletal muscle force model for the prediction of both muscle forces and surrounding musculoskeletal tissue loading for reasonably simplified anatomical models. While the skeletal muscle force model makes a significant contribution towards remedying the issue of indeterminacy in joint systems, it is unable to completely address the issue due to the number of forces that currently cannot be measured in vivo including joint contact loading, ligament loading, and tendon loading. For this reason, the anatomical model needs careful consideration in its development along with the exercise protocol for initial calibration of the muscle constant values, Km. These constants are a direct result of the dimensional analysis used to develop the model and represent a collection of the constant parameters that are influential in muscle force generation.1 The validation of the model’s solution process should also be considered under more complex conditions. For example, when observing the knee joint, the increased number of muscle forces and simultaneous use of both Newton–Euler equations (1) and (2) results in a much larger system of equations when utilizing several data points from each exercise. This increases the risk of finding local minima in the solution set which may not represent physiologically relevant solutions. Therefore, an analysis on the knee complex serves both as a validation of the model’s performance on a more complex joint as well as a validation on the solution process.
This study applied the skeletal muscle force model developed by Knodel et al.1 to the knee joint in order to characterize its performance on a more complex joint system where accuracy, and specifically subject-specific accuracy, is important. The solution process was also evaluated to determine its capability of delivering physiologically-relevant results. The final objective was to observe the calculated musculoskeletal tissue loading assuming that the first two questions were answered satisfactorily.
2. Methods
2.1. System setup and application of constitutive equation
The two Newton–Euler equations of rigid body motion were used in this analysis for observing the muscle forces and the surrounding musculoskeletal tissue loading and are given by

Fig. 1. Free Body Diagram of the Knee Joint Complex. The êi coordinate system’s origin is located at the joint center, J, and is fixed to the tibia. →FGM, →FGL, →FTFL, →FST, and →FBF are the muscle forces crossing the cut (red). →FLCL, →FMCL, →FACL, →FPCL (excluded from analysis), and →FPT are the tendon and ligament forces crossing the cut (blue). →FJM, →FJL and →FJ,LAT are the joint contact loads (yellow) where →FJM and →FJL represent the axial loads due to the medial and lateral condyles, respectively, and →FJ,LAT is the lateral joint contact load on the joint center. →FJ,LAT is originally assumed to be acting in the lateral direction for the purpose of writing out the equations. If the load comes out negative in the solution process, it is made positive and called →FJ,MED, acting in the opposite direction as shown in the FBD. Dashed portions of the force vectors represent sections of the force which are located behind the body of the tibia. Note that the full FBD encompasses the foot as well, meaning the ground reaction force, →FGRF, is acting on the foot outside of the scope of this figure. Similarly, the center of mass of the tibia, where gravity acts on the system, cannot be seen on this figure.
To reduce the number of unknowns, the muscle Km values introduced by Knodel et al. were assumed to be constant over time with the exception of fatigue and long-term physiological changes.1 Therefore, as more data points were observed from each exercise, the number of unknowns associated with the muscle force variables did not increase. All the other tissue loadings were able to vary with time. By removing the time dependency of the muscle forces from the equations, this application of the constitutive law effectively allowed a single one of the six scalar equations to cover the muscle forces given the amount of data points in time that were observed, leaving room for five time-dependent unknowns to be included in the anatomical model.
The five muscoloskeletal tissues included in the model were the anterior cruciate ligament (ACL), two axial contact loads due to the medial and lateral condyles, a medial/lateral joint contact load at the joint center, and either the medial collateral ligament (MCL) or lateral collateral ligament (LCL) (Fig. 1) based on the motion of the subject at the current point in time being analyzed. For instance, if the kinematic terms from (2) less the moment due to the ground reaction force resulted in a net positive moment about the ê1 direction, it was assumed the MCL was currently active. If the net moment was negative, the LCL was assumed to be active. This approach was taken because the MCL and LCL were not expected to both be carrying significant loads simultaneously. This would only happen if the tibia were being forcibly pulled away from the femur or if significant torque were occurring about the ê3 axis, neither of which occurred in the exercises from this study. Furthermore, assuming only one collateral ligament was carrying a significant load provided a reasonable simplification to the model when only five unknowns could be included outside of the muscle forces. The PCL was excluded from this analysis as posterior translation of the tibia was not expected in the conducted exercises. Additionally, PCL injuries can occur due to external contact to the knee,12 which were not experienced during these exercises, while ACL injuries are typically caused by non-contact actions where the tibia rotates laterally during knee flexion.13 Understanding the muscle force profiles will help elucidate the forces within the knee joint for a variety of scenarios, including non-contact ACL injuries on artificial turf.14 Such scenarios often combine frontal loading and rotation about the long axis of the tibia. The ACL also has a higher injury rate12 and so identifying its loading for these non contact exercises was more desirable. The geometry of the knee joint necessitated two unique axial joint contact loads for a more accurate loading distribution. A medial/lateral joint contact load was chosen for the last free unknown because it was expected to be more significant than an anterior/posterior joint contact load, especially for the cossack squat exercises. Relative translations at the tibiofemoral joint were left out of this analysis under the assumption that these minor translations would not have a significant impact on the muscles. This will, however, be an interesting point of investigation in future studies geared towards identifying ACL and PCL strains.
Five of the eight measured muscles directly cross the cut between the tibia and femur and two of those five are treated as known, specifically the two gastrocnemius muscles. This analysis is an extension of the validation study for this model done on the ankle joint2 which already calculated the KGM and KGL values of the subjects in this study. The three measured quadriceps muscles cross the cut by means of the patellar tendon. Due to the difficulty in tracking the kinematic information of the patella, which is needed to observe and solve for a true FBD of the patella, it was assumed that the magnitude of the patellar tendon load was equivalent to the magnitude of the combined quadriceps muscles. Lastly, only the inertia of the tibia was included in this analysis as the inertia of the foot was assumed to be negligible in comparison to both the kinematics of the tibia and the kinetics of the entire system.
2.2. Data collection
Five male and five female subjects of recreational athletic acumen between the ages of 21 and 26 participated in this study (Table 1). Purdue University’s Institutional Revfiew Board (IRB) approved the study and each participant signed a written consent prior to collecting data. Motion capture data was collected at 100Hz from eight Vicon Verio 2.2 cameras. An instrumented Bertec treadmill collected ground reaction force data at 1000Hz and Delsys Trigno EMG sensors sampled voltage readings of muscle activation at 2000Hz. The following muscles were measured with EMG for this analysis: the gastrocnemius medialis (GM), gastrocnemius lateralis (GL), rectus femoris (RF), vastus lateralis (VL), vastus medialis (VM), biceps femoris (BF), semitendinosus (ST), and tensor fasciae latae (TFL). Two-dimensional (2D) proton density (PD) magnetic resonance images (MRIs) were also taken of each subject’s right knee using a GE Discovery MR750 3T machine.
Subject | Age | Sex | Height [cm] | Mass [kg] |
---|---|---|---|---|
1 | 21 | Male | 177.8 | 63.5 |
2 | 26 | Female | 175.3 | 84.0 |
3 | 22 | Male | 179.1 | 64.9 |
4 | 22 | Male | 182.9 | 72.6 |
5 | 21 | Female | 165.1 | 61.2 |
6 | 21 | Male | 167.6 | 56.7 |
7 | 22 | Female | 170.2 | 63.5 |
8 | 20 | Female | 160.0 | 61.2 |
9 | 21 | Male | 177.8 | 74.8 |
10 | 26 | Female | 177.8 | 72.6 |
Optimal surface EMG (sEMG) placements were made by testing three nearby, separate locations on each muscle belly and choosing the location with the strongest signal. Maximum EMG signals were obtained prior to conducting the main exercise routine in order to acquire measurements for the Vref term in the constitutive equation.1,2 These readings were collected by having each subject perform restrained, explosive exercises for each muscle group. Resistance bands were used for leg extensions, hamstring curls, and hip abductions in order to target the quadricpes, hamstrings, and TFL, respectively. A weighted vest was used for single-legged calf raises for targeting the GM and GL (10–15% body weight based on subject ability).
Each subject performed ten repetitions of six unique exercises: a weighted and unweighted original squat (OS-W and OS-U, respectively), a weighted and unweighted cossack squat (CS-W and CS-U, respectively), and a weighted and unweighted lunge (L-W and L-U, respectively). Weighted versions of the squats and cossack squats were done with a weighted vest (20–25% body weight for squats and 15–20% body weight for cossack squats based on subject ability). Dumbbells were used for the weighted version of the lunges (15–20% body weight based on subject ability). Weighted and unweighted versions of each exercise were used in order to provide a wider range of muscle activation so that the Km values were better estimated.
2.3. Anatomical model
The constitutive equation from Knodel et al.1 requires various anatomical kinematic parameters including muscle moment arms and unit vectors. In order to estimate these parameters through time for each exercise, OpenSim3 was used. Each subject was first scaled in OpenSim and then inverse kinematics were performed on each exercise to generate motion files. These motion files were read in by a custom-built MATLAB script which then extracted the muscle attachment points, muscle lengths, joint angles, and the positional data of the tibia’s center of mass. The muscle attachment points were used to calculate muscle moment arms and unit vectors whereas the center of mass was used in estimating the inertia tensor of the tibia.
The MRI scans of each subject were used to obtain the patellar tendon’s moment arm and pennation angle as well as the unit vectors and moment arms of the ACL and PCL. While the PCL was excluded from this particular model, its kinematic data was still recorded for potential future use. These data points were extracted by inspecting the scans in an open source medical image viewer and using the built-in dimensional tools for identifying relevant distances between points of interest.
2.4. Estimating constant parameters
Various constant parameters are present in the constitutive equation which must be estimated prior to solving for each muscle’s Km value. This process has been described before in detail.1,2 In general, the parameters were originally estimated by either observing previously reported data and their corresponding curves or by performing nonlinear regression on the curves. They were then optimized by placing bounds around these estimates and identifying the values corresponding to the least amount of error.
2.5. Analysis
All data was collected and analyzed in a MATLAB program (Fig. 2) that evaluated the Km value of each muscle as well as the muscle forces, joint contact loading, and relevant ligament loading at each exercise peak. The details of this analysis are largely identical to that of the model’s validation study done on the ankle joint and can therefore be found in the corresponding paper.2 There are, however, some differences and the first is the identification of whether the MCL or LCL is engaged for the point in time being observed. This influences the population of the matrices as each equation is dependent on which of the two ligaments is currently carrying a load. The population of the matrices also differs between the two studies in that all six equations from (1) and (2) are being utilized simultaneously whereas the analysis of the ankle joint only initially used (2). This is due to the two axial joint contact loads not acting through the joint center, thereby generating moments about it and contributing towards (2).

Fig. 2. Analysis Flowchart derived from Fig. 5 of Knodel et al. (Ref. 2). *ϵ, c1, αasy, β2, a3, b3, β3.
Next, a simplified model was used to estimate the ratio between the quadriceps muscle force and the hamstrings muscle force so that bounds could be placed to ensure physiological relevance. Due to the complexity of the system and the necessary simplifications made, there are times when the solution process may find a local minimum corresponding to a zero muscle force. However, it is known that in many cases this is not reasonable and so the simplified model was used to place constraints. The simplified model consisted of the patellar tendon, a combined hamstrings muscle whose unit vector and moment arm were averaged from the two individual hamstring muscles (BF and ST), the two axial joint contact loads, the MCL or LCL depending on activation, and the medial/lateral joint contact load. The two gastrocnemius muscles (GM and GL) were also included since their Km values were already known from the previous ankle analysis. This simplified model was solved using the original squat exercises as they had the simplest motion and cleanest raw motion capture data. The ratio of quadriceps to hamstrings force was then calculated for each repetition and averaged. Constraints were placed on the solution process in MATLAB to ensure that the total quadriceps muscles and total hamstrings muscles maintained this ratio within +∕− 30%. It is not expected that this further simplified model predicts highly accurate quadriceps or hamstrings forces, but it is expected to predict a reasonable ratio between the two such that when the full model’s solution process is constrained to a set of bounds surrounding this ratio (e.g. +∕− 30%), it will produce physiologically relevant forces.
A second constraint was enforced on the quadriceps muscles due to their relation to the FBD via the patellar tendon. The earlier described relation of the four quadriceps muscle force magnitudes equating to the patellar tendon’s loading magnitude created an issue in the solution process. All four quadriceps muscles were dependent on the same unit vector and moment arm of the patellar tendon, rendering the program unable to identify the unique distribution of force among the four muscles. While the muscle model partially helps with this because of its dependency on each muscle’s current length, contraction velocity, etc., it does not completely remedy the issue. For this reason, it was assumed that the three measured quadriceps muscle forces were within +∕−30% of each other to avoid a single muscle dominating the other three, which was not probable for the exercise set used in this study. To estimate a force for the vastus intermedius (VI), which was not measurable via sEMG, a small portion of each of the other three quadriceps muscle forces were redistributed to the VI such that the patellar tendon force calculated by the mechanics remained the same. The remainder of the analysis was the same as the ankle analysis2 except that the joint contact loading was calculated alongside the Km values. This is a byproduct of using both (1) and (2) simultaneously whereas in the ankle analysis the joint contact loading was calculated after the Km values by observing (1). The OS-W exercise was chosen as the testing set for this analysis.
3. Results
The minimum and maximum Km values across the ten subjects in units of (C2⋅skg⋅m2∗109), respectively, were: 0.799 and 208.5 for the RF, 0.091 and 61.85 for the VL, 0.812 and 21.77 for the VM, 17.0 and 871.9 for the BF, 9.636 and 1786 for the ST, and 0.0027 an 52.72 for the TFL (Table 2). The significant variation in these values was also seen in a similar analysis at the ankle joint and is not surprising given the Km values’ initial dependency on EMG placement, which is accounted for by the model in subsequent sessions.2 Therefore, it is difficult to draw definitive interpretations from variations in a population’s Km values. Instead, they should be made from the calculated muscle forces.
Subject | KRF | KVL | KVM | KBF | KST | KTFL |
---|---|---|---|---|---|---|
1 | 2.403 | 0.752 | 5.324 | 140.6 | 68.73 | 5.042 |
2 | 14.96 | 4.648 | 4.599 | 871.9 | 1786 | 4.520 |
3 | 19.17 | 11.57 | 0.905 | 225.2 | 439.6 | 0.065 |
4 | 16.66 | 1.703 | 1.810 | 148.5 | 617.5 | 0.068 |
5 | 0.799 | 12.62 | 0.912 | 113.1 | 240.1 | 1.703 |
6 | 11.51 | 0.795 | 0.812 | 186.4 | 9.636 | 0.0027 |
7 | 6.513 | 0.091 | 1.597 | 17.00 | 18.06 | 0.0045 |
8 | 37.98 | 25.04 | 13.04 | 357.5 | 114.0 | 28.67 |
9 | 208.5 | 60.13 | 3.711 | 631.5 | 1634 | 49.62 |
10 | 54.49 | 61.85 | 21.77 | 800.8 | 792.2 | 52.72 |
Average percent errors were calculated in the same manner as the ankle analysis2 in that the kinematic components from (1) and (2) less the components of the forces and moments due to the ground reaction force and mass of the tibia were compared to the components generated by the predicted Km values (Figs. 3(a) and 3(b)). Average percent errors of all these components were calculated for each subject and ranged from 9.4% to 26.5% with an overall average of 20.6% across all subjects (Table 3).
Subject | Average percent error (%) |
---|---|
1 | 24.4 |
2 | 18.5 |
3 | 22.7 |
4 | 9.4 |
5 | 25.8 |
6 | 16.5 |
7 | 21.4 |
8 | 24.8 |
9 | 15.8 |
10 | 26.5 |
Average | 20.6 |

Fig. 3. Subject 1 Comparison of Predicted and Measured →B Components. Each set of bars represents a specific combination of direction (ê1, ê2, ê3) and repetition (i.e. ê1 of repetition 1 followed by ê2 of repetition 1 and so on, left to right).
Joint contact loading and ligament loading was also calculated by means of (1) and (2). When the lateral joint contact load was calculated as negative, it was renamed as FJ,MED and made positive. Average loads were calculated for each exercise and loads that were dependent on activation (i.e. MCL vs LCL and medial vs lateral JCL) had averages calculated only on repetitions in which the load was active (Figs. 4(a), 4(b), 5(a) and 5(b)). The MCL typically had larger loads than the LCL across all exercises and subjects and had higher average exercise loads than the LCL for each of the six exercises. The ACL had more variability in its loading between exercises and subjects and had an overall average loading that was comparable to the MCL. The axial JCLs dominated the medial/lateral JCLs and had similar magnitudes between each other across all exercises excluding OS-W. They also had the most variability between subjects (Table 4).

Fig. 4. Ligament loading as ratio of percent body weight. Recorded loadings occurred at the end of the descent phase for each exercise.

Fig. 5. Joint contact loading as ratio of percent body weight. Recorded loadings occurred at the end of the descent phase for each exercise.
Exercise | Average | Standard | ||||||
---|---|---|---|---|---|---|---|---|
Loading | CS-U (%) | CS-W (%) | L-U (%) | L-W (%) | OS-U (%) | OS-W (%) | (%) | deviation (%) |
ACL | 17.7 | 50.8 | 11.6 | 49.7 | 34.7 | 61.6 | 37.7 | 40.5 |
MCL | 27.9 | 29.0 | 29.7 | 31.2 | 41.9 | 40.3 | 33.3 | 14.8 |
LCL | 15.3 | 12.6 | 13.2 | 12.1 | 7.8 | 11.5 | 12.1 | 7.2 |
FJM | 175.9 | 238.9 | 145.3 | 223.9 | 124.3 | 184.8 | 182.2 | 58.0 |
FJL | 156.3 | 273.0 | 121.0 | 244.4 | 155.9 | 264.6 | 202.5 | 68.0 |
FJ,LAT | 0.9 | 5.2 | 1.5 | 4.9 | 7.3 | 9.4 | 4.9 | 6.6 |
FJ,MED | 11.2 | 10.7 | 9.2 | 9.4 | 2.2 | 2.3 | 7.5 | 6.8 |
Muscle forces were calculated for each subject for each point in time that was observed in the analysis, specifically the peaks of each repetition. Average muscle forces were calculated for each exercise across these peaks and represented as percentages of subject body weight (Figs. 6 and 7). Muscle forces were also calculated over time for each repetition and an averaged repetition was found for each exercise (Fig. 8). The quadriceps muscle forces were typically the largest in each of the exercises with the hamstring muscles accounting for either the second largest magnitude or a magnitude similar to that of the quadriceps muscles. The TFL and calf muscle forces were typically negligible in comparison to those of the quadriceps and hamstrings at the peak of the repetitions except in a few isolated incidents. One such example was in subject 1 who had considerable GL forces for the cossack squat and lunge exercises (Figs. 8 and A.1).

Fig. 6. Subject 9 muscle forces as percentages of body weight. Calf muscles are represented in shades of green, quadriceps muscles are represented in shades of blue, and hamstring muscles are represented in shades of orange. *The VI was estimated based on the total patellar tendon force as it was not measurable via sEMG.

Fig. 7. Averaged muscle forces as percentages of body weight across all subjects. Calf muscles are represented in shades of green, quadriceps muscles are represented in shades of blue, and hamstring muscles are represented in shades of orange. *The VI was estimated based on the total patellar tendon force as it was not measurable via sEMG.

Fig. 8. Subject 3 muscle forces over averaged exercise repetitions. For the lunge exercises, repetition starts at heel strike.
4. Discussion
The goals of this study were to further evaluate the predictive ability of the model by Knodel et al.,1 specifically at a more complex joint system, validate the solution process by observing its ability to produce physiologically relevant results, and identify the various surrounding tissue loadings at the knee joint during the exercises conducted for this analysis. The average percent error across all subjects was 20.6% (Table 3) and muscle forces were considered reasonable given the exercises and points in time where the muscle forces were observed (Figs. 6 and 9). Surrounding musculoskeletal tissue loading had the most variability in the ACL and two axial JCLs. The MCL, on average, carried a higher load than the LCL when engaged and was similar in magnitude to the ACL. In future studies, the model will be further refined by leveraging motion capture to investigate knee adduction and abduction in order to identify which ligaments are currently engaged rather than using resultant kinetics to drive the anatomical model. The two axial JCLs averaged approximately two times body weight across all exercises, with the CS-W generating the largest load for each. The medial/lateral JCLs, in comparison, were substantially smaller (Table 4).

Fig. 9. Subject 1 Semitendinosus muscle force vs. EMG signal for single repetition of L-U exercise. This plot demonstrates the change in both muscle force and EMG signal over the course of a repetition and that the force at the end of the repetition is not always the peak force for a given muscle or exericse.
The 20.6% average error is 12.4% higher than that of the ankle complex analysis.2 However, it was expected that the analysis at the knee would produce higher errors due to the anatomical complexities and the number of load-bearing structures necessitating simplifications of the anatomical model. For instance, only two of the four ligaments could be active at a time in this analysis due to the restriction on the number of unknowns that could be included in the system. Similarly, only three JCLs were permitted in the system and the geometries of the femur and tibia enforced two separate axial loads. A single axial load acting through the joint center was deemed to be too great of a deviation from the true loading distribution and would have compromised the model. Given this context, only one further JCL could be present and, for the reasons described earlier in the paper, a medial/lateral load was chosen over a posterior/anterior load. While all of these simplifications to the model were necessary and resulted in a plausible model, it was nevertheless missing the full description of the joint system. The ankle analysis resulted in simplifications to the anatomical model as well; however, they resulted in more mild alterations and therefore were expected to have less intrinsic error associated with the model. Additionally, the majority of the error in the →B components corresponded to the torque about the ê3 axis (Fig. 3(b)), which was considerably smaller than those about the primary ê2 axis (flexion/extension). In fact, if the torque components about ê3 were ignored for subject 1, the average percent error in the testing set dropped from 24.4% to 3.0%. When applied to all subjects, the total average percent error drops from 20.6% to 8.8%. This closely resembles the ankle analysis average of 8.2%2 and seems to indicate that the limitations of the anatomical model are the main contributor towards the error in the torque in the ê3 direction rather than the skeletal muscle force model, which performed similarly on the knee complex as it did on the ankle complex. It also demonstrates the importance in the development of an appropriate anatomical model for each study and its corresponding exercises, motions, or movements. In fact, the incorporation of MRI-generated geometry for the anatomical model should provide the most appropriate model on a subject-by-subject basis, when available. Results also indicate that the solution process was able to converge to a physiologically relevant solution as the quadriceps muscle forces were typically the largest and demonstrated co-contraction with the hamstrings muscles (Fig. 6). This phenomenon is expected as the hamstrings are responsible for stabilizing the joint during the descent of each repetition in all of these exercises. Occasionally, the hamstrings and quadriceps would share similar loads (Fig. A.1, OS-U), typically for the cossack squat or lunge exercises. The contributions of the calf muscles and TFL were negligible in most cases (counter-example, see Fig. A.1) as expected considering there was no plantar flexion or hip abduction for any of the exercises performed. The medial/lateral JCLs were shown to be much smaller in magnitude than the other musculoskeletal tissue loadings (Table 4). However, the loading in their relevant exercises indicate that they remain significant to the model. For example, the medial JCL was consistently near 10% body weight for the cossack squats and lunges and the lateral JCL was 7.3–9.4% body weight for the original squats while the remaining medial/lateral JCLs were substantially smaller. The MCL loading was consistently larger than the LCL loading for all exercises. This was not surprising as it suggests that the ground reaction force was typically located medially of the joint center, which is expected for each of the six exercises. However, there were still cases where notable LCL activation occurred (Table 4), meaning that either the ground reaction force was acting near the joint center and occasionally reversed between acting medially and laterally of the joint center, or there was enough error in the scrubbing and pre-processing of the data to affect which collateral ligament was considered active. The two axial JCLs increased during the weighted version of each exercise compared to the unweighted version, as expected. This was also true for the ACL but not for the MCL or LCL, suggesting that the ACL may be dependent on applied load while the collateral ligaments are not. This may be reasonable as the collateral ligaments are more prone to varus and valgus rotations, which would not necessarily increase due to an increased load. Instead, they would be dependent on the motion of the subject. Therefore, unless the increased load was large enough to alter the exercise motion, the MCL and LCL would not be expected to increase with exercise loading. The other potential reasoning for this behavior would be that the added loads were simply not large enough to warrant increased MCL and LCL activation.
The tensile strengths of the MCL and LCL are approximately 800N and 400N, respectively, and their yield strengths are 655N and 361N, respectively.15 The tensile strength of the ACL in cadavers has been reported across a wide range from 600–2300N13,16,17,18,19,20 while it is estimated that the tensile strength of the ACL in healthy adults is approximately 2000N.21 These tissue properties are encouraging as they provide upper extremes of what to expect from the calculated ligament loadings of this study. After converting the ligament loadings from percent body weights (Table 4) to pure loads (Table 5), it can be seen that all of the ligament loadings fall well below their respective tensile strengths. Furthermore, the MCL and LCL loading also fall well below their approximate respective yield strengths. All of this is expected as these exercises should not have produced maximum, or near maximum, loading conditions in the ligaments.
Exercise | ||||||||
---|---|---|---|---|---|---|---|---|
Tissue loading | CS-U | CS-W | L-U | L-W | OS-U | OS-W | Average [N] | Tensile strength [N] |
ACL | 109.9 | 308.2 | 70.7 | 309.1 | 221.9 | 404.3 | 237.4 | 200021 |
MCL | 181.9 | 186.3 | 192.3 | 205.1 | 277.6 | 261.6 | 217.5 | 80015 |
LCL | 99.1 | 83.7 | 85.5 | 77.2 | 48.3 | 74.2 | 78.0 | 40015 |
The primary limitation of this analysis is the simplified version of the anatomical model, as described earlier. These simplifications included the removal of any loads due to the PCL or an anterior/posterior JCL as well as the restriction of simultaneous MCL and LCL activation. The ACL was chosen over the PCL because of its higher rate of injury from non-contact situations13,12 and because posterior translation of the tibia was not expected during the selected exercises. The geometry of the knee articulation suggested that two separate axial loads be present in the model due to the medial and lateral condyles. A single axial load through the joint center would not only misrepresent the load distribution on the femur and tibia, but would also be unable of capturing the torque created by the axial loads about the ê1 axis. This torque is not insignificant given the magnitude of the calculated axial JCLs (Table 4). Lastly, a single collateral ligament was chosen to be active at any given time as neither a separating force between the femur and tibia nor a significant rotation of the tibia relative to the femur was anticipated. For these reasons, the model presented is expected to provide a reasonable simplification of the naturally complex knee joint and yield overall accurate results for the exercises conducted. The results were validated by the errors which suggest that any additional error present in this analysis relative to the ankle analysis are found in the comparatively small torques of the ê3 component. Another limitation arises due to the simplified anatomical model, more specifically the model of the patella. The lack of kinematic data of the patella during the exercises made it impossible to analyze the patella using mechanics. The aforementioned solution of essentially treating the patella as a pulley system with the patellar tendon and quadriceps muscles made finding the proper distribution of force between the four quadriceps muscles more difficult. For this reason, while total quadriceps force is expected to be accurate, the distribution of the total force among the RF, VL, VM, and VI has some uncertainty associated with it. If kinematic data of the patella can be obtained reasonably well and without too much difficulty for future studies, this would create more confidence in the distribution of force as a separate FBD of the patella could be incorporated into the overall model. Other limitations are similar to those of the ankle analysis2 including error accumulation and a generalized anatomical model in terms of bone geometry and corresponding muscle attachment points. Error is accumulated throughout the entire process of this analysis starting with the data collection (primarily with the motion capture data when markers were dropped) and followed by the scrubbing and pre-processing of the data prior to calculating the Km values. Furthermore, any errors that were present in the Km values of the ankle analysis were carried through to this analysis. This suggests that it would be difficult, if not impossible, to achieve better predictive ability of the model on a subsequent joint system compared to the original joint system. The anatomical model in this analysis was also generalized in terms of the bone geometry and muscle attachment points. While subject-specific patellar tendon and ACL kinematics were obtained, OpenSim was used to scale a generalized lower extremity model for observing muscle attachment points through time. However, it is known that there are significant variations in biological systems,22 so this simplification will lead to errors in the location of muscle attachment points and the subsequent calculation of muscle moment arms. Lastly, the proposed muscle model does not yet consider tendon dynamics.23,24 After examining the sensitivity of the current model and identifying feasible simplifications, future work will involve incorporating tendon dynamics, including both length and deformation velocity, into the model by considering these physical parameters in the initial dimensional analysis. Once these parameters have been incorporated, the model will be compared to other popular muscle modeling techniques, such as forward dynamics, to gain additional insights to the modeling process.
The results of this study have demonstrated the predictive ability of the skeletal muscle force model presented by Knodel et al.1 on the more complex knee joint system. They also validated the solution process by producing physiologically relevant muscle forces for the exercises performed in the study. Tracking patella kinematics could further enhance the solution process’ ability to properly distribute the quadriceps force among the four corresponding muscles. Furthermore, the joint system’s surrounding musculoskeletal tissue loading appears to behave as expected and is reasonable in magnitude. These findings are indicative of the muscle force model’s ability to predict muscle forces and tissue loading at any joint system assuming that a reasonable anatomical model can be created and the primary muscles are measurable via sEMG. For increased performance quality, the incorporation of MRIs into the anatomical model for the purpose of extracting subject-specific muscle attachment points should result in stronger predictive ability. This model could be used for both biomechanical investigations as well as for subject-specific studies where accurate feedback is important for injury prevention and performance optimization.
ORCID
ERIC A. NAUMAN https://orcid.org/0000-0002-8907-5716
Acknowledgments
We thank Norvin Bruns of Purdue University for assisting with the design and construction of the Vicon Vero camera mounting system. We also thank Dr. Bradley Duerstock of Purdue University for loaning his Delsys Trigno EMG system to make this study possible.
Appendix

Fig. A.1. Subject 1 muscle forces as percentages of body weight. Calf muscles are represented in shades of green, quadriceps muscles are represented in shades of blue, and hamstring muscles are represented in shades of orange. *The VI was estimated based on the total patellar tendon force as it was not measurable via sEMG.