|Home | About | Journals | Submit | Contact Us | Français|
Methods to predict contact stresses in the hip can provide an improved understanding of load distribution in the normal and pathologic joint. The objectives of this study were to develop and validate a three-dimensional finite element (FE) model for predicting cartilage contact stresses in the human hip using subject-specific geometry from computed tomography image data, and to assess the sensitivity of model predictions to boundary conditions, cartilage geometry, and cartilage material properties. Loads based on in vivo data were applied to a cadaveric hip joint to simulate walking, descending stairs and stair-climbing. Contact pressures and areas were measured using pressure sensitive film. CT image data were segmented and discretized into FE meshes of bone and cartilage. FE boundary and loading conditions mimicked the experimental testing. Fair to good qualitative correspondence was obtained between FE predictions and experimental measurements for simulated walking and descending stairs, while excellent agreement was obtained for stair-climbing. Experimental peak pressures, average pressures, and contact areas were 10.0 MPa (limit of film detection), 4.4-5.0 MPa and 321.9-425.1 mm2, respectively, while FE predicted peak pressures, average pressures and contact areas were 10.8-12.7 MPa, 5.1-6.2 MPa and 304.2-366.1 mm2, respectively. Misalignment errors, determined as the difference in root mean squared error before and after alignment of FE results, were less than 10%. Magnitude errors, determined as the residual error following alignment, were approximately 30% but decreased to 10-15% when the regions of highest pressure were compared. Alterations to the cartilage shear modulus, bulk modulus, or thickness resulted in ±25% change in peak pressures, while changes in average pressures and contact areas were minor (±10%). When the pelvis and proximal femur were represented as rigid, there were large changes, but the effect depended on the particular loading scenario. Overall, the subject-specific FE predictions compared favorably with pressure film measurements and were in good agreement with published experimental data. The validated modeling framework provides a foundation for development of patient-specific FE models to investigate the mechanics of normal and pathological hips.
It is estimated that 3 percent of all adults over the age of 30 in the United States have osteoarthritis (OA) of the hip , causing pain, loss of mobility and often leading to the need for total hip arthroplasty. Considerable clinical, epidemiological, and experimental evidence supports the concept that mechanical factors play a major role in the development and progression of OA [2-5]. For example, it has been demonstrated that a combination of duration and magnitude of contact pressures and shear stresses on the acetabular and femoral cartilage of hips with acetabular dysplasia can predict the onset of OA [6,7].
The ability to evaluate hip joint contact mechanics on a patient-specific basis could lead to improvements in the diagnosis and treatment of hip OA. To this end, both experimental and computational approaches have been developed to measure and predict hip contact mechanics (e.g. [6,8-14]). Experimental studies have been based on either in vitro loading of cadaveric specimens [8,10,13,14] or in vivo loading using instrumented femoral prostheses implanted in live patients [11,15-17]. While in vitro experimental studies have provided baseline values of hip joint contact pressure, testing protocols are inherently invasive, mechanical data are limited to the measurement area and specific joint pathologies cannot be readily studied. The use of instrumented prostheses represents the current state of the art for experimental study of in vivo hip mechanics [11,15-17]. However, the method is highly invasive and existing data are from older patients who have already been treated for advanced OA. There are no experimental methods available to assess hip contact mechanics non-invasively on a patient-specific basis.
Computational modeling is an attractive alternative to experimental testing since it is currently the only method that has the potential to predict joint contact mechanics non-invasively. Prior computational approaches have included the discrete element analysis (DEA) technique [18-20] and the finite element (FE) method [9,12,21]. These models have proven useful in the context of parametric, phenomenological or population (generalized patient model) based investigations. However, their ability to accurately predict patient-specific contact mechanics is questionable due to over-simplification of joint geometry and an absence of model validation [18-20].
Before computational models can be applied to the study of patient-specific hip joint contact mechanics, it is necessary to demonstrate that the chosen modeling strategy can produce models with accurate predictions and that the sensitivity of model predictions to variations in known and unknown model inputs is quantified . To fulfill these requirements we developed a subject-specific hip joint FE model (i.e. an FE model that is created from a single cadaver hip joint). Specifically, the objectives of this study were: 1) to develop and validate a subject-specific finite element (FE) model of hip joint contact mechanics using experimental measurements of cartilage contact pressure under physiological loading, and 2) to assess model sensitivity to several measured and assumed model inputs.
A combined experimental and computational protocol was used to develop and validate a subject-specific FE model of a 25 y/o male cadaveric hip joint (body weight = 82 kg). The joint was screened for osteoarthritis and the cartilage was determined to be in excellent condition (Tonnis Grade 0) .
All soft tissue with the exception of articular cartilage was removed. The acetabular labrum was dissected free from cartilage. Kinematic blocks were attached to the femur and pelvis for spatial registration between FE and experimental coordinate systems . The blocks were used to define anatomical axes for referencing joint loading angles using Bergmann’s coordinate system definition [15,16]. A volumetric CT scan of the hip was obtained (512 × 512 acquisition matrix, 320 mm field of view, in plane resolution = 0.625 × 0.625 mm, 0.6 mm slice thickness) using a Marconi MX8000 CT scanner (Phillips Medical Systems, Bothell, WA). The femur was dislocated from the acetabulum to ensure separation between the acetabular and femoral cartilage in the image data. A solid bone mineral density phantom (BMD-UHA, Kyoto Kagaku Co., Kyoto, Japan) was included to correlate CT voxel intensities with calcium equivalent bone density [25,26]. The aforementioned scanner settings produce thickness errors of less than 10% for simulated cartilage  and bone  when geometry is at least 1.0 and 0.75 mm thick, respectively.
Experimental loading was based on published data for in vivo hip loads [15,16]. Bergmann et al. reported hip joint anatomical orientations (flexion, abduction, rotation) and equivalent hip joint forces (magnitude and direction) during routine daily activities for 4 patients with instrumented femoral prostheses [15,16]. Data for their “average patient” (mean data for the 4 patients analyzed) were used in the present study to simulate walking, stair-climbing and descending stairs. A custom loading apparatus was developed to apply the kinematics that corresponded to these loading conditions (Fig. 1).
The iliac crests of the pelvis were cemented into a mounting pan in neutral anatomical orientation (anterior iliac spine in plane with plane of pubis symphysis, [15,16]) and attached to a lockable rotation frame (Fig. 1 left). The rotation frame was flexed and abducted relative to the vertical axis of the actuator to simulate the orientation of the equivalent hip joint force vector for each loading scenario. The femur was potted and attached to a lockable ball joint (Fig. 1 middle). Three-dimensional orientation of the joint was achieved by flexing, abducting and rotating the femur relative to the pelvis. Equivalent joint reaction force angle and anatomical orientation were confirmed by digitizing the loading fixture surfaces (joint force) and planes of the kinematic blocks (anatomical orientation) using a Microscribe-3DX digitizer (Immersion Corp., San Jose, CA) with a measured positional accuracy of ±0.085 mm . The digitized points were fit to planes, and angles between the planes were calculated. The orientation of the pelvis and femur fixtures were adjusted until the direction of the joint reaction force vector and anatomical orientation angles were within ± 1° of those reported by Bergmann’s average patient.
Low range (1.7-10 MPa)-pressure sensitive film (Sensor Products Inc., Madison, NJ) was used to measure joint contact pressures. Preliminary testing assessed the efficacy of two other ranges of film before choosing the low range film. Both super-low (0.5-1.7 Mpa) and medium film (10+ Mpa) were tested. The super-low film was nearly completely saturated following loading, which indicated that nearly all of the recorded pressures were > 1.7 MPa. On the other hand, the medium range film (10+ Mpa) only indicated minuscule areas where pressures were > 10 MPa. Prior to dissection, different film sizes were cut into a rosette pattern using a knife plotter. The film size that maximized contact area and minimized overlap was chosen (Fig. 1 right). Small notches were cut in the anterior, posterior, and medial aspects of the rosette to reference the location of contact pressures relative to the hip joint.
Peak loads for each activity were simulated by displacing the femur into the acetabulum at a constant rate. For each activity, the rate of actuator displacement was adjusted until peak loads were achieved within 0.33 sec, representative of the time required by the average subject reported by Bergmann et al [15,16]. Three to four cycles of preconditioning were necessary to obtain the correct displacement rate. The pressure film was then attached to the head of femur between sheets of polyethylene. Planes of the kinematic blocks were digitized to establish an experimental coordinate system in neutral orientation. The femur was then displaced into the acetabulum until the target load was achieved. The actuator was returned to its starting position at the same displacement. The three notches on the film were digitized. The specimen was allowed to recover between trials for over 100 times the interval that was needed to reach peak load. The entire protocol was repeated 3 times for each loading scenario. The films were stored in a dark location for 48 hours following testing  and then scanned and converted to digital grayscale images. An independent calibration curve was established to relate pixel intensity to pressure .
Commercial software was used to segment surfaces of the cortical bone, trabecular bone, cartilage and kinematic blocks in the CT image data (Amira 4.1, Mercury Computer Systems, Boston, MA). Splines representing the outer surface of cortical bone were obtained from automatically thresholded images . Cartilage was segmented from air using a threshold value that resulted in the greatest accuracy for reconstructing simulated cartilage in a phantom based imaging study . The boundary between trabecular and cortical bone was segmented both automatically and semi-automatically. When cortical and trabecular bone blended together in the image data they were manually separated.
Triangular surfaces were generated for each structure from the segmented image data using the Marching Cubes algorithm . The outer cortical surface facets were decimated to achieve a polygon surface discretization that was consistent with our previous study . Cartilage surfaces were decimated and smoothed slightly to remove visible triangular irregularities and segmentation artifact. The triangular surface mesh for cortical bone was converted to a quadratic 3-node shell element mesh [25,32-34]. Position dependent shell thickness was assigned to each node, based on the distance between adjacent trabecular bone boundary nodes . The resulting pelvis and femur cortical meshes consisted of 13,562, and 4,196 elements, respectively (Fig. 2), representative of the mesh density that has been shown to produce accurate predictions of cortical bone strains in prior pelvic FE modeling . The interiors of the cortical shell meshes were filled with tetrahedral elements to represent trabecular bone of the femur and pelvis . The final pelvis and femur trabecular bone tetrahedral meshes consisted of 227,108 and 82,176 elements, respectively, which was consistent with our prior study of pelvic FE modeling .
Acetabular and femoral cartilage surfaces were imported into FE preprocessing software (TrueGrid, XYZ Scientific, Livermore, CA) and hexahedral element meshes were created. Convergence studies were performed by increasing the number of elements through the thickness of the cartilage incrementally while the overall aspect ratios were held constant by adjusting the in-plane mesh resolution. The meshes for acetabular and femoral cartilage were considered converged if there was less than a 5% change in contact area, peak pressure and average pressure between subsequent meshes.
Cartilage was represented as an incompressible, neo-Hookean hyperelastic material  with shear modulus G=6.8 MPa . Incompressibility was enforced using the augmented Lagrangian method . Cortical bone was represented as hypoelastic, homogenous, and isotropic with elastic modulus E=17 GPa and Poisson’s ratio ν=0.29 . Trabecular bone was represented as isotropic hypoelastic with ν=0.20 . An average elastic modulus was calculated for each tetrahedral element using empirical relationships from the literature [25,26] and the BoneMat software . Overlap between the shell and tetrahedral elements  was accounted for by assigning an elastic modulus of 0 MPa to all tetrahedral elements that shared nodes with shell elements.
To establish the neutral kinematic position for each loading scenario, the FE model was transformed from the CT coordinate system to the appropriate experimental reference system . Nodes superior to the pelvis cement line, those residing within the sacroiliac (SI) and pubis joint, and those inferior to the cement line of the potted femur were defined rigid according to anatomical boundaries determined experimentally. The rigid femur nodes were constrained to move only in the direction of applied load, while the nodes at the pubis, SI and cement line were fully constrained. The Mortar method was used to tie acetabular and femoral cartilage to the acetabulum and femoral head, respectively [40,41]. Contact between the femoral and acetabular cartilage was enforced using the penalty method . All analyses were performed using NIKE3D .
Sensitivity studies were performed to investigate how changes in assumed cartilage material properties, thickness and FE model boundary conditions affected predictions of cartilage contact mechanics. The baseline cartilage shear modulus was altered by ± 1 SD (G=10.45 and 2.68 MPa) using standard deviations for human cartilage . To ascertain the effects of the assumption of cartilage incompressibility, bulk to shear modulus ratios of 100:1 (ν=0.495) and 10:1 (ν=0.452) were analyzed. To account for differences in segmentation threshold intensity between real and simulated cartilage , the baseline threshold value used to segment cartilage was adjusted by ± 50%. Updated cartilage FE hexahedral meshes were generated based on these surfaces. To quantify the affects of model boundary conditions, three separate cases were analyzed: 1) bones were assumed rigid, 2) the rigid constraint at the pubis joint was removed, and 3) trabecular bone was removed so that deformation of only the cortical bone was considered. Separate models were generated for each loading activity, yielding a total of 27 models.
A program was developed to compare FE predicted cartilage contact pressures with results from pressure sensitive film. The program allowed for the investigation of two types of error: 1) misalignment between FE and experimental results, and 2) differences in the magnitude of contact. First, the program converted the grayscale images of pressure to fringed color using the calibration curve. Next, FE pressure predictions were transformed into a synthetic film image with the same dimensions, including rosette cuts, as the pressure films. Surface nodes of the femur cartilage FE mesh were fit to a sphere and then flattened by a spherical-to-rectilinear coordinate transformation. The synthetic image was aligned with the pressure film image using the experimentally digitized notches on the pressure film. The rosette cuts on the experimental film were duplicated in the synthetic FE pressure film image by moving the FE pressure results circumferentially, according to the wedge angle of the rosette. Separate synthetic FE images were created and aligned for each experimental image since the pressure films were not placed in the exact same anatomical position between loading trials. Finally, a difference image between each synthetic and experimental image was created by subtracting the corresponding pixel intensities.
A root mean squared (RMS) error criterion was used to assess the degree of similarity by comparing pixel intensity values between FE and experimental images. Only those pixel intensities within a user specified range were compared. This range was taken as the full sensing range of the film (1.7-10 MPa) but was also determined for smaller 2 MPa bins of pressure to assess the ability of the FE models to predict pressures within specific ranges. Further constraints were made in the calculation of RMS error because, in this study, the experimental film data were considered the “truth”. Specifically, if a pixel in the synthetic FE image indicated a pressure within the specified range but the corresponding experimental film pixel did not, then the pixel was not included in the calculation. However, if an experimental pixel was within the range but its corresponding FE pixel was not, then the pixel was included.
Misalignment error was also distinguished from magnitude error. Misalignment error could occur due to fitting the FE mesh to a sphere, from inaccurate digitization of the notches used to align the results, or from FE model inaccuracies. Misalignment error was quantified independently by calculating the RMS error real time while the synthetic FE image was rigidly rotated about a spherical coordinate system. The rotations required to minimize RMS error, along with the re-calculated anterior, posterior, and lateral positions of the notches were recorded. Misalignment error was then expressed as the difference in pre and post alignment RMS error whereas magnitude error was taken as the post alignment error.
Peak pressure, average pressure, contact area, and center of pressure were calculated for each experimental film and synthetic FE image after the synthetic FE films were aligned with experimental images to minimize RMS error. FE peak pressure was determined by recording the maximum FE pressure value within the region of experimentally measured contact. Experimental peak pressures were calculated as the maximum experimental film pixel intensity. Pixel intensities that indicated pressures within the film range (1.7-10 MPa) were used to determine FE and experimental average pressures. Contact area was calculated by multiplying the number of pixels within the detectable pressure range of the film by the area of each pixel (0.0154 mm2). The center of pressure (COP) was found by determining the center of the image, which was weighted according to pixel intensity. The difference in the centers of pressure between images (FE synthetic film COP - experimental film COP) was expressed as the anatomical difference in anterior/posterior and medial/lateral positions over the detectable range of the film.
Similar analyses were performed for the sensitivity studies. First, a baseline image was created for each baseline FE model from which subsequent sensitivity study predictions were compared. RMS errors were calculated per above. Changes in peak and average pressure, and contact areas were calculated per the methodology discussed above, but a larger pressure range was used (0.5 - 10 MPa) since the pressure range was no longer limited by the film. For the cartilage sensitivity studies, percent changes in peak pressure, average pressure and contact area were reported as combined results from the three loading scenarios analyzed. For the boundary condition sensitivity studies, results were reported as percent changes with respect to each loading scenario.
Cortical bone thickness in the pelvis and femur averaged 1.8±0.8 and 2.9±2.3 mm, respectively. The resulting trabecular bone modulus in the pelvis and femur averaged 270±188 and 295±198 MPa, respectively. Three elements through the cartilage thickness were necessary to yield converged FE predictions. The final mesh for acetabular and femoral cartilage consisted of 15,000 and 23,415 elements, respectively (Fig. 2). Cartilage thickness in the acetabular and femoral cartilage meshes was 1.6±0.4 and 1.5±0.5 mm, respectively, as estimated using the cortical bone thickness algorithm  (Fig. 3). The final FE model consisted of ~400,000 elements (Fig. 2), and each analysis took on the order of 2 hours of wall clock time.
Experimental pressures ranged from 1.7-10.0 MPa (range of film detection). All pressure films recorded pressures at the upper limit of film detection (10 MPa). However, less than 5% of the total pixels fell into this category. FE predictions of peak pressure were 10.78, 12.73 and 11.61 MPa for walking, descending stairs and stair-climbing, respectively. Experimental average pressure and contact area ranged from 4.4-5.0 MPa and 321.9-425.1 mm2, respectively, while FE predicted average pressure and contact area ranged from 5.1-6.2 MPa and 304.2-366.1 mm2, respectively (Fig. 4).
The experimental pressure recordings revealed bi-centric patterns of contact during simulated walking and descending stairs and a more or less mono-centric pattern during simulated stair-climbing (Fig. 5 top row). Experimental pressure distributions were similar during simulated walking and descending stairs, with a horseshoe shaped bi-centric contact pattern directed anterorlaterally to posterormedially. When the femur was rotated internally during simulated stair-climbing the contact pattern was oriented in a lateral to medial direction (Fig. 5 top row). Overall, the magnitude and location of FE predicted contact pressures corresponded well with experimental measures. However, experimental bi-centric contact patterns during simulated walking and descending stairs were not predicted by the FE models (Fig. 5 middle row).
Patterns of FE-predicted contact varied with the loading activity (Fig. 6 bottom row). The majority of contact occurred along the lateral aspect of the acetabulum for all three loading activities (Fig. 6 bottom row). The contact area moved from anterior to posterior as the resultant load vector changed from shallow extension during descending stairs to more moderate flexion angles during walking and stair-climbing (Fig. 6 bottom row).
Difference images, calculated prior to manual alignment between FE synthetic and experimental films, further clarified the degree of qualitative agreement between FE synthetic and experimental films (Fig. 5 bottom row). Differences in contact pressure were greatest for descending stairs and were least during stair-climbing (Fig. 5 bottom row). Overall, misalignment errors were less than 7% (Table 1). The rotations and resulting translations of the experimental film fiducials required to minimize RMS errors were less than 3° and 5 mm for walking and stair-climbing but were substantially higher for the descending stairs case (22° and 9 mm) (Table 1).
Following manual alignment, residual RMS errors were on the order of 30%. As suggested by the difference images, errors were greatest for descending stairs and least for stair-climbing (Fig. 5 bottom row). When RMS error was plotted in 2 MPa pressure bins, RMS errors decreased to around 10-15% at the maximum bound of pressure analyzed (8-10 MPa). This finding indicates that FE models were best suited for predicting the higher stressed regions of cartilage, corroborating the good qualitative correspondence between FE synthetic and experimental films in these locations (Fig. 5).
Differences in center of pressure locations, as calculated over the entire film detection range, were less than 10 mm (Table 2). The smallest difference in the COP occurred for stair-climbing, while the largest difference occurred during descending stairs. In general, COPs for the FE models were directed more lateral (-) and anterior (+) to experimentally measured COPs.
Changes of ±50% to the shear modulus resulted in approximately a ±30% change in FE predictions of peak pressures, while changes in average pressure and contact area were around ±10% (Fig. 7 top left). Lowering the cartilage Poisson’s ratio from ν=0.5 to ν=0.495 did not have an appreciable effect (Fig. 7 middle left). However, a further decrease in the Poisson’s ratio to 0.452 resulted in a 25% decrease in peak pressures, while changes in average pressure and contact area were less than 10% (Fig. 7 middle left). Altering the thickness of femoral and acetabular cartilage (~ 10% change average cartilage thickness) resulted in less than a ±10% change in FE predictions (Fig. 7 bottom left). RMS differences between baseline and all cartilage sensitivity study results were approximately 6.5%, indicating that the spatial distributions and magnitudes of contact pressure did not change substantially.
Rigid bone models decreased computation times from ~2 hours to less than 10 minutes. Representing the bones as rigid structures affected both the magnitude (Fig. 7 top right) and spatial distribution (Fig. 8) of cartilage contact pressure, but the degree of effect depended on the loading activity. RMS differences in synthetic films between baseline and rigid bone models averaged 29.2±5.5%. FE predictions of peak pressure, average pressure and contact area were altered but also varied according to the loading scenario analyzed (Fig 7 top right). When the rigid constraint on the pubis joint was removed, FE predictions changed on the order of −15 to +5% (Fig. 7 middle right). Finally, when the trabecular bone was removed, i.e. only the cortical shells supported the cartilage, changes in FE predictions ranged from −25 to +5% (Fig. 7 bottom right). Average RMS differences between baseline results and the latter boundary condition sensitivity studies were only 3.1%.
To our knowledge this is the first study to validate FE predictions of cartilage contact pressure with experimental measurements using a cadaveric human hip joint. The purpose of developing and validating a subject-specific model was to ensure that the chosen computational protocol could produce a model capable of predicting in vitro cartilage contact pressures. The FE model provided very reasonable predictions of both the spatial distribution and magnitude of cartilage contact pressure under the simulated loading conditions. Excellent predictions were obtained for simulated stair-climbing. The posterior aspect of the bi-centric experimental contact pattern was not predicted by the FE model for walking and descending stairs. Nevertheless, the magnitude of pressure in these locations was low in comparison to the anterior region where the FE models provided more reasonable correspondence.
Small manual rotations of the pressure film were necessary to minimize RMS errors for simulated walking and stair-climbing. In contrast, the descending stairs case required a substantial amount of manual rotation (Table 1). It is likely that the majority of misalignment error was due to the method of digitizing the film fiducials during the experiment. It was necessary to move the linear actuator up by ~20 mm to access the film markers. It was assumed that this displacement resulted in a perfect vertical translation for purposes of defining the marker coordinates, but when the coordinates were plotted relative to the translated model they did not reside on the surface of the cartilage. This offset was minor during walking and stair-climbing but was greater during descending stairs. The femur was in extension for this loading activity and when the translation was applied, the femoral neck would have contacted the edge of the acetabulum, resulting in an offset of the film marker coordinates. Contact in this location would not have occurred with the hip in moderate and deep flexion during walking and stair-climbing.
The finding that RMS magnitude errors decreased when the bounds of pressure were increased suggests that the models were best suited for predicting localized “hot spots”. Therefore, the modeling strategies developed herein may be well suited for predicting the primary region of contact, which may be sufficient for many patient-specific modeling applications.
FE predictions of average pressure and contact area were not overly sensitive to changes in the cartilage shear modulus, bulk modulus or thickness (±10%). However, greater changes in peak pressure were noted (up to ~25%). This finding demonstrates that peak pressure prediction requires more accurate model inputs for cartilage geometry and material properties than for average pressure prediction.
Computational models of the hip have often represented bones as rigid structures [12,44], which is an attractive simplification because solution times are greatly reduced. The present study demonstrated that the assumption of rigid bones can alter predictions of cartilage contact stresses in the hip. The effect is modulated by the specific boundary and loading conditions in the model. Because the consequence of the rigid bone assumption cannot be assessed without a direct comparison to the case of deformable bones, investigators should use caution when representing the bones as rigid for modeling cartilage contact mechanics in the human hip.
Although the contralateral pelvis was left intact in the experimental study, the FE models assumed that the pubis joint was rigid. The results of the sensitivity study that removed the pubis constraint demonstrated only minor differences in FE predicted cartilage contact mechanics, thereby giving credence to this model assumption. While this simplification was warranted for the current study, it may not be appropriate for models where load is directed more medially (e.g. simulations of side-impact loading ).
Since the reported elastic modulus of trabecular bone is orders of magnitude less than cortical bone, we investigated whether or not trabecular bone needed to be represented in the models. The results of the sensitivity study suggest that it plays little mechanical role with regard to cartilage contact stresses. Therefore, for patient-specific modeling applications it may be appropriate to exclude trabecular bone, assuming that similar boundary and loading conditions are assigned.
Experimental studies have used pressure sensitive film to measure hip joint contact pressures under similar loading conditions [8,13,14]. Peak pressures measured by von Eisenhart-Rothe et al.  ranged from 7 MPa at 50% body weight to 9 MPa at 300% body weight, in fair agreement with the results of the current study. Bi-centric, horseshoe shaped patterns extended from the anterior to posterior aspect of the femur were noted . Afoke et al.  measured peak pressures on the order of 10 MPa at 350% body weight and the anterorsuperior surface of the cartilage was identified as an area of high pressure . All of these studies reported irregular, non-symmetric pressure distributions [8,13,14].
Large differences in material properties, geometry and boundary conditions make it impossible to directly compare the FE predictions from this study with prior modeling studies, but some general trends can be identified. Nearly all FE hip joint modeling studies to date have used two-dimensional, plane strain models [9,12,21,44] with either rigid [12,44] or deformable bones [9,21]. To our knowledge, the earliest FE contact model was reported by Brown and DiGioia . In this study, FE predicted pressures were irregularly distributed over the surface of the femoral head. Values of peak pressure were on the order of 4 MPa at loads representative of those applied in the current study. Rapperport et al.  developed a similar model that predicted peak pressures on the order of 5 MPa at 1000 N of load. Using rigid bone models resulted in predictions only slightly different than the deformable bone model.
Macirowski et al.  used a combined experimental/analytical approach to model fluid flow and matrix stresses in a biphasic contact model of a cadaveric acetabulum. To our knowledge, this is the only previous FE study to explicitly model the acetabular cartilage thickness. The acetabulum was step loaded to 900 N using an instrumented femoral prosthesis, yielding peak contact pressures on the order of 5 MPa. When the experimentally measured total surface stress was applied to the FE model, average predicted pressures (solid stress + fluid pressures) were approximately 1.75 MPa. The lower range of pressure used to determine average pressures was not specified, making it impossible to directly compare average results. However, scaling the applied load of our model to 900 N and assuming a lower bound of 0.3 MPa to calculate average pressure (lowest pressure isobar indicated by Macirowski et al.) yields average pressures of 2.47±0.29 MPa over the three loading scenarios analyzed, which is in good agreement with Macirowski et al’s predictions.
Yoshida et al.  developed a dynamic DEA model to investigate the distribution of hip joint contact pressures using the Bergmann gait data. The model assumed spherical geometry and concentric articulation. Qualitatively, our predictions of primary contact during simulated walking, descending stairs, and stair-climbing are in good agreement with the results of this study, but the spatial distributions of contact were markedly different. Peak pressures during walking, descending stairs, and stair-climbing in the study by Yoshida et al. (3.26, 3.77, 5.71 MPa, respectively) were substantially less than those predicted in the current study.
With the exception of the study by Macirowski et al., the FE models developed herein predicted higher contact pressures than previous FE and DEA, studies. This discrepancy is most likely due to the assumptions of spherical geometry and concentric articulation in the prior computational studies. Although the literature suggests that normal hips may be modeled as spherical structures with concentric articulation [19,46], the hip joint is not spherical and cartilage thickness is not uniform [12,47,48].
The aforementioned computational models assumed a cartilage modulus ranging from 10-15 MPa [9,21] yet cartilage was given a baseline modulus of ~40 MPa (G=6.8 MPa) in the current study. While one might expect that a higher modulus would result in equivalently higher contact pressures, the results of our sensitivity studies demonstrate that this is not the case, as changes in the cartilage shear modulus of ±50% resulted in only a ±25% and ±10% change in peak and average contact pressures, respectively. Even with a 25% reduction, peak pressures predicted in this study were still nearly double those reported previously [9,18-20].
Several limitations of the current study must be mentioned. First, only one model was developed, based on a singe CT image dataset. However, extensive sensitivity studies were conducted on this single model to understand the importance of model inputs and material properties. In the future, we plan to develop several more subject-specific models to ensure that the protocol utilized herein produces accurate models. Secondly, experimental loads were based on average in vivo data from older patients who had already undergone treatment for advanced hip OA. Given the large inter-patient variation in joint kinematics observed by Bergmann et al. [15,16], the use of average loading data likely did not accurately represent the actual kinematics for the specimen in this study. Our approach is justified since the objective of the experimental protocol was to apply realistic loading and boundary conditions that could be reproduced in the FE simulations for model validation.
Pressure film was chosen to measure cartilage contact pressures in this study since: 1) it is reasonably accurate (10- 15% error [49,50]), 2) it can be cut into rosette patterns to conform to highly spherical surfaces (thereby preventing crinkle artifact) and, 3) it is inexpensive and has been used extensively in prior studies of hip joint contact stresses [8,13,14]. A limitation of film pressure measurement is that the technique records a “high watermark” rather than measurements of dynamic pressures or the time-loading history [51,52]. However, film measurements have been shown to be equivalent to the contact stresses resulting from an incompressible elastic analysis , making the use of pressure film appropriate in the current study. Pressure film would not be appropriate for dynamic loading studies (e.g. simulations over entire gait cycle) since it is inaccurate in shear. However, prediction of dynamic pressures was not a goal of this study. Rather, the objective was to predict cartilage pressures at peak joint reaction force in the gait cycle during static, unidirectional loading. Results from the pressure measurements indicate that contact occurred beyond the perimeter of the film during simulated walking and descending stairs. While it would be desirable to capture the entire region of contact, it was not feasible to do so using larger rosettes as they caused excessive overlap and crinkle artifact during pilot testing. Finally, it was found that all of the films were saturated making it impossible to determine the true values of experimental peak pressures. However, our pilot study demonstrated that the low range pressure film was the best choice for the experiment. In addition, less than 5% of the film pixels had saturated pressures, suggesting that peak pressures were very close to the saturation limit of the film (10 MPa).
Removal of the labrum is noted as a potential limitation to this study since this structure was healthy in the specimen tested. The labrum was removed because: 1) labral geometry could not be distinguished from the CT data as a separate entity since it had the same image intensity as the adjacent cartilage, and 2) labral tissue properties and a corresponding constitutive equation have not been extensively reported in the literature.
The literature regarding the contribution of the labrum to hip cartilage mechanics is unclear [44,54,55]. Using pressure film, Konrath et al.  found no significant changes in contact area, mean pressure, or maximum pressure in the anterior or superior acetabulum and only noted a significant decrease in the maximum pressure in the posterior acetabulum when the labrum was removed. In contrast, an in vitro study Ferguson et al.  demonstrated that hip joints with the labrum removed consolidated more and had substantially decreased intra-articular pressures under both constant and cyclical loads. As a precursor to including the labrum in future FE modeling studies, it is clear that more extensive material testing is necessary to characterize the labrum’s constitutive behavior along with effective methods to distinguish this structure from adjacent cartilage in the image data.
Although actions of individual muscles were not considered, the equivalent joint reaction force was based on in vivo data [15,16]. The primary focus of the present research was to quantify cartilage contact pressures in the peri-acetabular region rather than bone stresses in areas where muscles were attached. Therefore, we could justifiably model the action of all muscles as a single equivalent force vector acting through the hip joint.
Although cartilage exhibits biphasic material behavior , it was represented as incompressible hyperelastic in this study. In vitro studies suggest that fluid flow is minimal during fast loading [12,54], making our assumption of incompressibility warranted given the loading rates used in the experiments. We recently demonstrated the equivalence between biphasic and incompressible hyperelastic FE predictions during instantaneous loading . Cartilage also exhibits depth dependent material properties , variation in stiffness over its surface [43,59] and tension-compression nonlinearity . Incorporating these aspects might have resulted in different, perhaps better, predictions of contact stress magnitude and distribution. Future modeling efforts should investigate the importance of these effects via sensitivity studies.
As discussed above, simplified analytical models or population based approaches to patient hip joint modeling have not yielded predictions that are consistent with in vitro data. Although it is possible that these discrepancies are due to model assumptions, it is difficult to pinpoint sources of error unless some reference standard (i.e. experimental data) is available for comparison. Therefore, we believe subject-specific modeling is a necessary precursor to either population or patient-specific modeling. The benefit of using a subject-specific approach first is that computational predictions can be directly compared to data obtained experimentally. The ability to directly quantify model accuracy is lost in population or patient based studies as direct validation is impossible. With a validated protocol in place it becomes much more feasible to develop patient-specific models that provide clinically meaningful data in terms of improving the diagnosis and treatment of hip OA, and for the study of pathologies such as hip dysplasia.
In conclusion, our approach for subject-specific FE modeling of the hip joint produced very reasonable predictions of cartilage contact pressures and areas when compared directly to pressure film measurements. Predictions were in good agreement with other experimental studies that used pressure film, piezoelectric sensors and instrumented prostheses [8,12-14,61,62]. The sensitivity studies established the modeling inputs and assumptions that are important for predicting contact pressures. The validated FE modeling procedures developed in this study provide the basis for the future analysis of patient-specific FE models of hip cartilage mechanics.
Financial support from NIH #1R01AR053344 and F31EB005551 is gratefully acknowledged.