Patellofemoral pain (PFP), a frequent complaint in orthopedic practice, is associated with 25%–40% of knee injuries (Liao et al., 2015). This pain is aggravated by various activities that load the joint (e.g., squatting, running, ascending/descending stairs), greatly limiting the daily activities of individuals with PFP (Fick et al., 2022). There is no consensus on the exact mechanism of PFP development (Vannatta and Kernozek, 2015). However, a prevailing theory holds that PFP develops in response to increased patellofemoral joint stress (PFJS) (Salsich and Perman, 2007). Chronic overuse of the patellofemoral joint has been related to increased intraosseous pressures and pain, microfractures, increased bone metabolism, and increased bone water content, all of which have pathologically detrimental effects on the subchondral bone (Ho et al., 2014). Studies have reported that individuals with PFP exhibit elevated PFJS during walking and running when compared with those in pain-free individuals (Farrokhi et al., 2011b; Liao et al., 2015). Based on the associations, PFJS seems to be an important factor in assessing the load on the patellofemoral joint. It can be helpful in preventing patellofemoral joint injuries, evaluating the effectiveness of PFP rehabilitation programs.

The analytical model is currently the most widely used method for evaluating the PFJS, which is based on the formula obtained from previous cadaver experiments. The classic analytical model has several significant limitations, including the failure to account for the synergistic contraction of knee joint muscles and the consideration of only sagittal plane factors in relation to PFJS (Bonacci et al., 2014; Atkins et al., 2019). To achieve greater accuracy in assessing PFJS, various methods have been developed, including musculoskeletal models, discrete element analysis (DEA), and finite element analysis (FEA). Nunes et al. systematically reviewed the literature which was utilized analytical models to evaluate PFJS and described the possible best paradigm to evaluate PFJS (Nunes et al., 2018). An important limitation of the review is its relatively narrow definition of the methods of PFJS. Nunes et al. do not take into consideration the methods of musculoskeletal modeling, DEA, and FEA. It is necessary to carry out a new literature review due to the progress of PFJS evaluation technology and the increase in the number of studies since the publication of the previous review.

The primary goal of the present systematic review was to identify and categorize the methods developed and used to evaluate PFJS comprehensively by taking into consideration all important aspects of PFJS (e.g., synergistic contraction, variables of the coronal plane, and 3D geometry of bones). This paper will help researchers fully understand the assessment schemes for patellofemoral joint stress and different optimization strategies, while also proposing future research directions.

While conducting this review, we followed the recommendations of the PRISMA statement (Supplementary Table S6) (Liberati et al., 2009).

Studies of the patellofemoral joint mechanism have particularly concerned the forces and stresses in the joint (van Eijden et al., 1986). Non-invasive evaluation of PFJS will provide important information on the factors contributing to joint load in a special population and may provide data that can be used to guide treatment (Brechter and Powers, 2002b). Thus, it is necessary to evaluate the PFJF or PFJS non-invasively. This review aims to identify and categorize methods used to evaluate PFJS non-invasively, and we found four methodologies: analytical model musculoskeletal model, DEA, and FEA.

At present, there is still no “gold standard” for PFJS assessment because it is impossible to measure the PFJS directly in vivo. Accurate assessment of PFJS will aid engineers in the design of better tissue-engineered constructs for cartilage replacement and assist physiotherapists to evaluate the PFP patient’s recovery and design a more effective intervention. Researchers have adopted some techniques to optimize the PFJS assessment program to improve the accuracy of PFJS or expand its scope of application. It is mainly aimed at patella kinematics and quadriceps muscle strength optimization.

The geometry or alignment of the patella affects kinematics, contact mechanics, and strain in the patellar bone, potentially leading to PFP (Brechter and Powers, 2002b; Fitzpatrick and Rullkoetter, 2012). The differences in subject-specific anatomy among the population will cause errors in the process of calculating the PFJS. Researchers use personalized programs to eliminate errors. In the analytical model evaluation scheme, the patellofemoral contact area was estimated using a function of knee flexion angles from 0° to 90°, as previously described in past cadaver studies. It is not clear whether the cadaveric data can represent the characteristics of the patellofemoral joint in vivo. Furthermore, in studies of patellofemoral bracing, it is not appropriate to calculate the patellar joint contact area by formula because the patellofemoral bracing will affect the patellar trajectory and change the patellofemoral contact area (Powers et al., 2004a). This implies that the generalized equation for the patellofemoral contact area may not be fully accurate (Zavala et al., 2021). A analytical model based on MRI technology was applied to establish a regression equation of contact area based on their respective sample (Brechter and Powers, 2002a; b; Powers et al., 2004a; Whyte et al., 2010). However, the MRI assessment protocol still has some limitations. For example, participants completed the scan in a relaxed quadriceps state to avoid motion artifacts in the MRI, which may not reflect the patellofemoral joint contact area when the quadriceps is contracted.

The accurate prediction of patellofemoral joint kinematics is a key prerequisite prior to investigating the effects of evaluating the impact of injury and physical activity on the PFJS. The geometry of the patella and femur is an important factor affecting the kinematics of the patella. Yu et al. (2019) study revealed that soft tissue plays a critical role in adjusting patellar tracking during the early stages of knee flexion, while the morphology of the trochlear groove and patellar face determines the relative position of the patella and femur as the knee flexion angle increases. Among the eight included studies that used musculoskeletal models to evaluate PFJF or PFJS, two musculoskeletal models were used. The first is the standard model available in simulation software (e.g., AMMR 1.6.2, MoCapModel) (Alexander and Schwameder, 2016), including patella, femur, and tendon, and the model was scaled to match each participant’s anthropometry and mass (Andersen et al., 2010). The second is a subject-specific musculoskeletal model based on subject-specific bone geometry, cartilage surface, and muscle fiber orientation (from MRI) (Chen and Powers, 2014; Lenhart et al., 2015a). Similarly, the included studies of PFJS evaluation using DEA and FEA also adopted the subject-specific model based on the radiation image. Using standardized models to evaluate PFJS may lack biofidelity. The subject-specific model may contribute to more accurately predicting patella tracking and PFJS evaluation (Mesfar and Shirazi-Adl, 2005).

During active knee flexion, the force of the quadriceps muscle is transmitted to the patellar tendon to generate stress on the contact surface of the patellofemoral joint. The researchers used several methods to estimate quadriceps muscle strength. The simplest method is to divide the knee extension moment by the effective lever arm of the quadriceps muscle (from radiography or magnetic resonance imaging). This method does not consider the co-activation of the hamstring and quadriceps during the activity. The second method used net joint moments and muscle moment arms (function of joint angles) to derive hamstring, quadriceps, and gastrocnemius muscle forces (Willson et al., 2015a; Willson et al., 2015b; Willy et al., 2016; Sinclair et al., 2018; Starbuck et al., 2021; Zavala et al., 2021). The third method is to modify the quadriceps muscle forces by the muscle’s cross-sectional area, maximum voluntary contraction force per unit cross-sectional area, and the EMG-force relationship (Escamilla et al., 2008a; Escamilla et al., 2008b; Escamilla et al., 2009). The fourth method is to estimate quadriceps muscle strength by summing the knee flexion moment and the net moment calculated by inverse kinematics. The knee flexion moment was estimated by SIMM software based on the individual’s lower extremity kinematics, the velocity of muscle contraction, and flexor muscle EMG (Chinkulprasert et al., 2011; Powers et al., 2014; Teng and Powers, 2014). The fifth method is computed muscle control algorithm that modifies muscle excitations by using feedforward and feedback control to follow recorded joint angle trajectories (Thelen et al., 2014). In addition, the researchers used static optimization techniques to estimate quadriceps muscle strength. Specifically, an active force adjustment, passive force adjustment, and velocity adjustment were interpolated from Arnold et al. (2010) curves for a given fiber length and velocity. The maximal dynamic muscle force was then calculated from the maximal isometric muscle force (Kernozek et al., 2015). From the perspective of the calculation principle, the quadriceps muscle force calculated by the coordinated contraction adjustment, computed muscle control algorithm, and static optimization techniques, while not flawless in their ability to represent muscle force in vivo, is more accurate than the force calculated from the net moment of knee extension and moment arm. However, it is still difficult to accurately solve the muscle force because the musculoskeletal system is redundant (the identical joint moment can frequently result from an unlimited number of combinations of muscle forces) (Johnson et al., 2022). And the motion simulation depends on experimental data and various uncertainty parameters, such as the variability of marker position, motion artifact, how to normalize EMG, and collecting EMG from deep muscles (Ackermann and van den Bogert, 2010; Burden, 2010). An earlier study has shown that the muscle is activated but does not produce effective muscle force in individuals with PFP (Herzog, 2000). Muscle force is a major factor affecting PFJS, and inaccurate quantification of muscle force may lead to misinterpretation of contact stress patterns. Islam et al. used three-dimensional registration techniques and linear mapping to measure the virtual penetration depth from the patellar cartilage surface to the femoral cartilage surface, combined with a FEA to estimate PFJS (Islam et al., 2015). The method provided by Islam et al. does not require the calculation of muscle force but provides accurate kinematics data of the patella and femur. However, the scanning space of MRI technology is narrow, and it is difficult to capture dynamic functional activities (walking, running, etc.). The development of biplane fluoroscopy systems has provided a method for directly measuring the kinematic data of the patella relative to the femur during functional activity (Wheatley et al., 2020). Gustafson et al. (2021) used a DEA model driven by high-precision kinematic data collected by a biplane fluoroscopy system to evaluate the PFJS during walk tasks. The method of DEA driven by high-precision kinematic data can estimate the PFJS at a subject-specific level without estimating the muscle force, avoiding the source of error. Therefore, the method for evaluating PFJS has a theoretically high accuracy. At the same time, this method can provide information on the change in the stress curve and distribution of the patellofemoral joint during the entire functional activity period, thereby promoting the understanding of PFP and the innovation of rehabilitation methods. The methods described are applicable for studying individuals with abnormal patellar trajectories and additional knee joint load.

FEA and DEA are refined computational models that allow the integration of subject-specific musculoskeletal parameters and experimental data, and are of great value for understanding stress distributions in complex biological structures (Fernandez et al., 2008; Farrokhi et al., 2011b; Gustafson et al., 2019). The researchers believe that shear stress is related to cartilage degeneration and damage, and higher shear stress in cartilage may also be the cause of damage to the exciting nociceptors present in the subchondral bone plate of the patella (Wojtys et al., 1990; (Besier et al., 2008). The FEA has been used to calculate the shear stress of PFP patients and healthy people during activities, and it was found that PFP patients have a higher level of shear stress (Farrokhi et al., 2011b). The studies of FEA and DEA included in this review conducted MRI or CT scans to obtain specific bone geometries to build subject-specific models. The material properties of cartilage are determined using constitutive models of varying complexity. These models range from linear spring models, which are described by a single constant, to isotropic, linear elastic models with two material constants, and a biphasic model with three material constants (Keenan et al., 2013). The studies used a linear elastic material model, which is widely accepted as valid due to the cartilage’s short-term elastic response during activities with loading frequencies over 0.1 Hz, such as walking (Higginson and Snaith, 1979). The average element size in the FEA was 0.5–1 mm, and previous studies on convergence analysis have shown that element sizes smaller than 0.75 mm did not result in significant stress changes but resulted in longer simulation times (Farrokhi et al., 2011b; Liao et al., 2015). The cartilage material properties of Young’s modulus and Poisson’s ratio are important factors affecting stress distribution of cartilage. Unfortunately, researchers do not seem to have reached a consensus about the material properties of cartilage. The minimum elastic modulus used in the included literature is 2 MPa, and the maximum is 40 MPa. Indentation testing is the main means of measuring the properties of materials. Although the indentation test is very accurate for determining material properties, researchers have not yet reached a consensus on cartilage properties. The thickness of cartilage, dehydration of the cartilage, contact area of the indentation test, and strain rate will affect the stiffness of cartilage. Kempson et al. reviewed the study of measuring the cartilage elastic modulus and gave the range of the cartilage elastic modulus under different working conditions (Kempson, 1980). Although the elastic modulus used in most studies related to PFJS evaluation is within this range, the elastic modulus used in some studies is beyond this range, and the basis for the elastic modulus was not provided. It should also be noted that changes in cartilage thickness will affect biomechanical properties (Shaktivesh et al., 2019). Individuals with PFP have been confirmed to have a possible loss of cartilage thickness (Farrokhi et al., 2011a). The rationality of using the same elastic modulus as that of healthy people to evaluate the PFJS of individuals with PFP needs further study. However, there is currently no research directly assessing the mechanical properties of articular cartilage in patients with PFP. The main technical limitation is the inability to directly estimate the mechanical performance of articular cartilage through non-invasive assessments such as MRI (Brenneman Wilson et al., 2023). Indirect prediction methods based on cartilage composition, such as quantitative MRI (Hatcher et al., 2017) and Raman spectroscopy (Bonifacio et al., 2010) and numerical simulations for performance prediction (Gupta et al., 2009), have shown promising potential in non-invasively assessing the mechanical properties of cartilage. In future research, these indirect prediction methods can be used to predict the mechanical properties of cartilage in patients with PFP, thereby promoting the development of stress calculation models for the patellofemoral joint.

In summary, the evaluation of PFJS using analytical models remains the most commonly used approach. Currently, researchers are enhancing the accuracy of PFJS evaluations by combining analytical models with MRI or musculoskeletal models. Musculoskeletal simulations, FEA, and DEA offer more detailed information on stress distribution. FEA can compute stress variations throughout the depth of the cartilage and estimate shear stress, which may be a key indicator of pain. Compared to analytical models, these three approaches consider more PFJS influencing factors and reduce sources of error. However, they are more time-consuming and require technical expertise. The evaluation of PFJS is undoubtedly a complex task. If the goal is to achieve high accuracy in evaluating PFJS, combining advanced imaging technologies (such as biplane fluoroscopy system) with advanced computational techniques (DEA and FEA) would be a good option. The method should be further developed to establish the “gold standard” for PFJS evaluation. If the focus is on exploring the impact of changes in muscle strength on PFJS, utilizing musculoskeletal models would be a suitable approach. For those seeking a simpler computational solution, analytical models would be appropriate. However, there is no clear boundary in the application scope of the aforementioned approaches. In practical applications, the appropriate technical means and evaluation approaches need to be selected based on theoretical considerations and/or ecological validity. Furthermore, the insufficient consideration of confounding factors remains a key limitation to the current research quality. There are multiple factors that impact PFJS such as sub-group classification of PFP patients (Selfe et al., 2013) and foot posture (Neal et al., 2014). In future studies, it is necessary to exercise better control over variables and minimize the influence of confounding factors. These measures will contribute to an improved level of evidence in research.

This study has potential limitations. A total of 69 studies were included in this review and the PFJF and PFJS data of included studies were extracted, but it was not possible to quantitatively analyze the calculation results between different methods due to the high heterogeneity (population, task).

There are four methods to evaluate patellofemoral joint reaction force and patellofemoral joint stress, including the analytical model, musculoskeletal model, discrete element analysis, and finite element analysis. At present, there is still no “gold standard” for PFJS. And researchers are still trying to improve the evaluation accuracy of PFJS. This is mainly achieved through using a personalized model and optimizing (or avoiding) quadriceps muscle strength calculation. In theory, the evaluation scheme of combining advanced computational and biplane fluoroscopy techniques has high accuracy in evaluating PFJS. Additionally, in practical applications, it is important to select the appropriate technical methods and evaluation approaches based on theoretical considerations and/or ecological validity. In the future, researchers can develop relevant model frameworks to accurately calculate PFJS and provide technical solutions for a better understanding of the mechanism of patellofemoral joint pain and optimization of the patellofemoral joint treatment program.