Evaluation and Validation of Thorax Model Responses: A Hierarchical Approach to Achieve High Biofidelity for Thoracic Musculoskeletal System

As one of the most frequently occurring injuries, thoracic trauma is a significant public health burden occurring in road traffic crashes, sports accidents, and military events. The biomechanics of the human thorax under impact loading can be investigated by computational finite element (FE) models, which are capable of predicting complex thoracic responses and injury outcomes quantitatively. One of the key challenges for developing a biofidelic FE model involves model evaluation and validation. In this work, the biofidelity of a mid-sized male thorax model has been evaluated and enhanced by a multi-level, hierarchical strategy of validation, focusing on injury characteristics, and model improvement of the thoracic musculoskeletal system. At the component level, the biomechanical responses of several major thoracic load-bearing structures were validated against different relevant experimental cases in the literature, including the thoracic intervertebral joints, costovertebral joints, clavicle, sternum, and costal cartilages. As an example, the thoracic spine was improved by accurate representation of the components, material properties, and ligament failure features at tissue level then validated based on the quasi-static response at the segment level, flexion bending response at the functional spinal unit level, and extension angle of the whole thoracic spine. At ribcage and full thorax levels, the thorax model with validated bony components was evaluated by a series of experimental testing cases. The validation responses were rated above 0.76, as assessed by the CORA evaluation system, indicating the model exhibited overall good biofidelity. At both component and full thorax levels, the model showed good computational stability, and reasonable agreement with the experimental data both qualitatively and quantitatively. It is expected that our validated thorax model can predict thorax behavior with high biofidelity to assess injury risk and investigate injury mechanisms of the thoracic musculoskeletal system in various impact scenarios. The relevant validation cases established in this study shall be directly used for future evaluation of other thorax models, and the validation approach and process presented here may provide an insightful framework toward multi-level validating of human body models.


INTRODUCTION
The human thorax is a fairly complex body region with musculoskeletal structures (ribcage) and soft tissue organs in the thoracic cavity. The ribcage consists of 12 thoracic vertebrae (T1-T12), 12 pairs of ribs, costal cartilages, sternum, scapulae, and clavicles. In automobile crashes, thoracic injury is one of the most prevalent injuries, and ranks second only to head injury with reference to the number of fatalities and serious injury outcomes (Forman et al., 2019). According to the Abbreviated Injury Scale (AIS), an anatomically based injury severity scoring system, the musculoskeletal (MSK) injury accounted for the highest percentage of thoracic AIS in all passenger vehicle occupants (Cavanaugh and Yoganandan, 2015). The thoracic MSK system protects intrathoracic structures from external penetration or blunt trauma, which is the leading cause of injury during motor vehicle collisions (MVC). Biomechanically, the external force under blunt impact is mostly sustained by the bony skeleton of the thorax, and the impact energy can be dissipated by muscles and soft organs with viscoelasticity properties. A better understanding of the thorax biomechanics, especially the biomechanical response and injury mechanisms of thoracic MSK, is crucial to prevention and mitigation of thorax injuries.
To investigate MVC related human injuries and enhance vehicle safety designs, surrogates of occupants have been used, including post mortem human surrogates (PMHSs), anthropomorphic test devices (ATDs), and volunteers. They each have their strengths and limitations, but their common feature is that they provide a similar response to that of the living human body during a crash and are studied to assess MVC injuries and develop or improve design for automotive safety systems (Beeman et al., 2012). Over the past three decades, computational finite element (FE) human body models (HBMs) have emerged as a powerful and versatile tool that can accurately describe the human anatomy, biomechanics, and variability for injury risk predictions and vehicle safety systems development. Compared to surrogates for safety designs, computational thorax FE models can integrate and represent thoracic biomechanical data from PMHS or volunteers, as well as efficiently evaluate both overall thoracic responses and localized physical variables related to injuries. For example, FE models can quantitatively predict complex chest deformation and bone fracture patterns that are not able to be evaluated via ATDs (Pipkorn and Kent, 2011).
A considerable number of thorax FE models have been developed to simulate chest behaviors under complex kinematics during various impact loading conditions. Yang et al. (2006) reviewed thorax FE models published before 2005 and after 2005 (Yang, 2018). These models were either isolated thorax models (Deng et al., 1999;Shen et al., 2008) or integrated into full body model (FBM) (Robin, 2001;Haug et al., 2004;Kimpara et al., 2005;Song et al., 2009;El-Jawahri et al., 2010;Pipkorn and Kent, 2011;Antona-Makoshi et al., 2015;Iwamoto et al., 2015;Poulard et al., 2015;Schoell et al., 2015;Zhu et al., 2016). Despite the significant difference of modeling details, the primary goal of these models was to evaluate human thoracic response, understand thoracic injury tolerances, and assess injury risks under various impact loads. Ensuring the biofidelity of these FE models, particularly for the thoracic MSK system, is indispensable for the credibility and utility of the tool for predicting injury to the ribcage and the intrathoracic tissues and organs. In general, one of the major limitations of most of these thorax models was the lack of reasonable validations, mainly due to a significant paucity of available experimental data at multi-levels. For example, if a thorax model was validated under a specific impact condition, its kinematics and injuries might not be valid under different impact loads. Even if the model was globally validated, it is uncertain if it has the capability of predicting thoracic responses and injury risk at a localized level. Additionally, the tissue material parameters are not always available, and most of the experimental datasets at the component or organ level were not prepared for the purposes of numerical simulations and model validation. As such, HBMs have been developed by the Global Human Body Models Consortium (GHBMC), which is a multi-institution collaboration involving industry, government, and academic researchers aiming to develop the most biofidelic computational HBMs for crash safety advancement. Although the GHBMC models contain a significant amount of anatomical details and can provide a basis for further model development and applications, the thoracic region in the detailed GHBMC midsized male occupant FE model (M50-O) v4.5 model (released in 2018, the base model used for this study) has inadequacies and limitations, particularly in some important structures of the MSK system. These limitations could be attributed to model simplifications or a lack of material verification, which requires further enhancement to increase accuracy and model biofidelity.
The goal of this study was to report on the hierarchical improvement, evaluation, and validation of the GHBMC M50 thorax model, focusing on the biofidelity of the thoracic MSK system at multi-levels. Specifically, this study focuses on improvements and validations to the thoracic spine, clavicle, sternum, costal cartilage, etc. To our knowledge, FE models for thoracic spine (T-spine) are fewer compared to cervical and lumbar spines. Detailed T-spine models with proper verifications and validations (V&V) using experimental studies are relatively rare, especially for impact scenarios. Aira et al. (2019) developed and validated a simplistic T-spine model using quasi-static loading for adjacent vertebrae and rear hub impact for chest force-deflection response. However, it did not include the stateof-the-art modeling techniques and injury-predictive capabilities that have been used in C-spines models (Panzer et al., 2011). This study builds on these previous efforts to enhance the existing GHBMC T-spine model with improved biofidelity, additional validations, and injury simulation capability. The enhanced model was validated against several experimental cases at segment, functional spinal unit (FSU) and full torso levels. The T-spine was incorporated into the thorax model, and it was able to predict the soft tissue deformation and ligament injury response under traumatic loading levels. The costovertebral joints were evaluated at a component level across the whole ribcage. These substructures were validated for multiple directions of motion, including ventral-dorsal and cranial-caudal flexion. For bony skeleton tissues, the material parametric studies were implemented to determine the properties of the clavicles, sternum, and costal cartilages. Different failure thresholds using maximum principle strain (MPS) were added to these structures to characterize the bone failure features. Several experimental setups were accurately replicated to simulate the axial compression of clavicle, bending testing of clavicle and sternum, and cantilever-beam like bending testing of costal cartilages. The bone behavior of the ribs was not included here since they were already well investigated in previous studies Poulard et al., 2015;Zaseck et al., 2018).
To evaluate the biofidelity of the whole thorax model, body regional level and full thorax or torso level validation tests were conducted across different loading conditions. These tests contained point loading of the eviscerated ribcage, frontal pendulum impacts, shoulder pendulum impact and table-top belt loading tests. The biomechanical response of the thorax was overall deemed biofidelic and the model was found to be in good agreement with the experimental data using both qualitative and quantitative assessments. The validated thorax model, as demonstrated in this study, was able to serve as a valuable tool for safety researchers and automobile designers to predict, prevent, and mitigate thoracic injuries of vehicle occupants. Additionally, this study describes a multi-level, hierarchical evaluation and validation strategy to enhance the biofidelity of other human body models for injury risk prediction and prevention.

Body FE Model Overview
The GHBMC 50th percentile male (M50) detailed seated occupant (v4.5) FBM was the basis for the improvement, development, and validation in this study. The model was developed for use with the LS-DYNA solver, a commercial nonlinear explicit finite element analysis (FEA) program (LSTC, Livermore, CA, United States). The anatomical geometry of the body was reconstructed based on the anthropometric data scanned from a 50th percentile male volunteer (Gayzik et al., 2009). The high quality quadrilateral and hexahedral meshes of the thorax components were generated by an interactive multiblock meshing approach Poulard et al., 2015). The FBM contained 1,036 parts, which were discretized by 2.19 million of elements and 1.26 million nodes.
In the base thorax model, the bones of thoracic cage were modeled as elasto-plastic materials with strain-rate dependency ( * MAT_Piecewise_Linear_Plasticity) in LS-DYNA, including cortical and trabecular bones of the ribs, sternum, and clavicles. The thoracic vertebrae were represented by rigid bodies, and the costal cartilages were modeled as an elastic material. For each thoracic intervertebral disc, the annulus and nucleus pulposus was modeled as a highly compressible Hill foam (Hill, 1979) and an inviscid fluid, respectively. The thoracic spine ligaments were defined as non-linear elastic spring beams without failure properties. A supplemental table was provided to collect the material models and properties for the important parts relevant to the evaluation and improvement in current study (Supplementary Table 1).
The process of the hierarchical improvement, evaluation, and validation of the thorax model is presented in this section (section "Thoracic Intervertebral Joints Evaluation and Improvement, " section "Costovertebral Joints Evaluation and Enhancement, " section "Clavicle Model Improvement and Validation, " section "Sternum Model Improvement and Validation, " section "Costal Cartilage Model Enhancement and Validation, " section "Ribcage or Full Thorax Level Validation Cases, " and section "Model Performance Evaluation"), and a graphical overview of the steps of the study can be found in the Supplementary Figure 1).

Thoracic Intervertebral Joints Evaluation and Improvement
The T-spine was enhanced and improved by accurate representation of the components, material properties, and ligament failure properties at tissue level. For the intervertebral discs (IVDs), the nucleus pulposus and annulus fibrosus ground substance were modeled using solid elements. For each IVD, four pairs of concentric quadrilateral layers (8 layers total) were created to represent the annulus fibrosus fiber laminae (Newell et al., 2017), which were embedded in the ground substance ( Figure 1A; Panzer and Cronin, 2009). There were 9 major spinous ligaments of the T-spine included in the model, including anterior and posterior longitudinal ligaments, ligamentum flavum, interspinous ligaments, supraspinous ligaments, intertransverse ligaments (left and right), and facet capsular ligaments (left and right) ( Figure 1B). The ligament response was characterized using non-linear load-deflection curves identified by three distinct spinal regions (Chazal et al., 1985;Panzer et al., 2011). The failure deflection for each individual ligament was adopted from data reported in Pintar (1986).
Initial quasi-static evaluations of the model were conducted for adjacent vertebrae and the intervertebral disc between Frontiers in Bioengineering and Biotechnology | www.frontiersin.org them (motion segment) to ensure a valid tissue-level response (Markolf, 1972), including axial tension and compression, shearing and lateral bending. Since the experimental testing of motions segments was not able to represent natural spinal motions due to too many constraints imposed, the FSU with three adjacent vertebrae and two intervertebral discs was therefore adopted by researchers to evaluate more realistic spinal responses (Yang, 2018). In our study, two sets of FSUs (T2-T4 and T7-T9) were simulated (Figure 2A) and compared with experimental testing results of flexion bending (Lopez-Valdes et al., 2011). In addition to the FSU level evaluation, the full torso response and change of spine angle of extension were validated by tests of rear hub impact loading (Forman et al., 2015) with a 97.5 kg impactor using three different velocities: 1.5, 3.0, and 5.5 m/s, shown in Figure 2B.

Costovertebral Joints Evaluation and Enhancement
The costovertebral joints were evaluated at component level using the method outlined in Duprey et al. (2010). To replicate the experiment setup, the ribs 2, 4, 6, 8, and 10, along with their adjacent vertebrae were isolated from the thorax base model. All relevant ligaments that interconnect the bony structures were preserved for each level of the costovertebral joint, including the radiate, costotransverse ligament, superior costotransverse ligament, lateral costotransverse ligament, and the intra-articular ligament. The radiate was modeled using shells and assigned as a non-linear fabric material. Other ligaments were modeled as tension-only beams with a single constant to define the stiffness. Each rib was inserted by a simple model of a cylindrical rod, which was rigidly constrained to the rib. Four quasistatic rotational motions were applied to the rod to investigate ventral-dorsal and cranial-caudal flexion, as demonstrated in Supplementary Figure 2.

Clavicle Model Improvement and Validation
The clavicle model was evaluated at the component level using experimental testing results in Zhang et al. (2014), which loaded the clavicles to fracture by two methods: three-point bending and axial compression, respectively. The FE models were created based on an isolated clavicle from the thorax model and modeling for the testing device. For the three-point bending test, each end of the clavicle was supported by a pinned assembly which allowed rotation along the superior-inferior axis. The impactor in Figure 3A, which was modeled as cylinder with aluminum material properties, was loaded along the anteroposterior direction (i.e., perpendicular to the clavicle longitudinal axis) using the same rate as in the experiments (0.1 m/s). In the axial compression test, as shown in Figure 3B, the medial end was potted into a square-shaped block that allowed rotation, while the lateral end was clamped and moved along the lateral-medial direction by 0.1 m/s. In the experiments, four uniaxial strain gages were attached at the perimeter of the clavicle cross section, as indicated by green rectangles in Figures 3A,B. Based on the force-deflection responses, parametric studies of the material properties were conducted to determine the appropriate elastic modulus and yield strength of the cortical bone. Bone failure was defined using a maximum principal strain (MPS) criterion. Through a brute-force optimization process, it was found that 8% MPS for cancellous bone and 3% MPS for cortical bone produced reasonable results to the experimental data (Duprey et al., 2008). Once the MPS reached the defined threshold, the element was eroded and removed from the calculation.

Sternum Model Improvement and Validation
The sternum fracture was evaluated under a dynamic three-point bending test (Kerrigan et al., 2010). The test rig was modeled similar to the three-point bending test of the clavicle ( Figure 3C). The impactor was put in the middle of the two posts and displaced along the downward direction (1.115 m/s) to represent anteroposterior loading on the sternum. Parametric studies of material properties were conducted to determine appropriate sternum material parameters (elastic modulus, yield strength and failure strain) based on experimental force-time responses.

Costal Cartilage Model Enhancement and Validation
The costal cartilage model was evaluated at the component level using experimental testing data from Forman et al. (2010). Figure 4 showed the testing device with potted costal cartilage model, which was isolated from the segment connected to the anterior border of the forth rib. The end of cartilage connected to the rib (i.e., right side of the device) was constrained, and the left end was loaded along the vertical direction (0.4 m/s) to represent the situation of mid-chest compression. Anatomically, the costal cartilage is composited of an inner solid of hyaline cartilage and a surrounding tissue layer called the perichondrium. During FE modeling, the solid structure of costal cartilage was updated from a linear elastic material model in base model to an elasto-plastic material with strain-rate dependency. Additionally, a non-linear fabric shell layer that surrounded the solid mesh was incorporated to include the biomechanical contribution of the perichondrium, since it was not explicitly modeled previously. The material parameters were calibrated to the reaction force  along the anteroposterior direction between the simulation and the reported experimental data.

Ribcage and Full Thorax Level Validation Cases
After validating the MSK structures at the component level, the components were integrated into the thorax model in FBM. The overall mechanical response of the chest, which included the contribution of these MSK soft and hard tissues, could be evaluated against a number of loading cases at ribcage or full thorax level. Here four experimental cases were selected based on their relevance, the availability of data, and the large range of structural levels they represented.

Point Loading of the Eviscerated Ribcage
The quasi-static point loading of the ribcage was simulated using the experimental testing approach reported in Kindig et al. (2010). In the experiment, isolated eviscerated ribcages were mounted upright and quasi-statically loaded by a plate interfacing with a spherical segment glued to the superficial surface of the ribs or the sternum. The reaction force on this plate against applied displacement was measured to check the biomechanical response. To simulate the testing, the FE model of the ribcage (Figure 5) was positioned consistent with the initial position adopted in the experiments. The loading and boundary conditions were also set similarly to the experimental settings: the thoracic vertebrae were constrained and a constant velocity (0.2 m/s) was imposed on the plate up to the prescribed displacement used in the experiment (varied with the point loading site).

Frontal Pendulum Impact
The frontal pendulum impact testing was designed to quantify the thoracic response and injury tolerance along the anteroposterior direction under midsternal blunt impact. The PMHS used in the experiments was placed in a seated upright position, and compressed on the frontal chest by a 23.4 kg cylindrical impactor with an initial impact velocity of 4.3 m/s and a diameter of 152 mm centered with the sternum at the level of the 4th intercostal space. In the simulation, the FBM was seated on a rigid plate and the impactor positioned at the same mid-sternum level as in the PMHS testing (Figure 6). The overall thoracic responses (deflection-time and force-deflection relationships) were reported and compared with experimental corridors (Kroell et al., 1971;Lebarbe and Petit, 2012).

Shoulder Pendulum Impact
The lateral impact on the shoulder provided response characteristics and injury tolerance of the shoulder and upper thorax along the lateral direction. Koh et al. (2005) presented a meta-analysis of biomechanical responses from several cadaveric studies of lateral impacts on the shoulder, which were developed into biomechanical response corridors for four loading speed cases: 4.4 m/s (foam padded impactor). 4.5 m/s (unpadded), 6.4 m/s (padded), and 6.8 m/s (unpadded). To replicate the experiments, the model was seated on a rigid plate with the arms down by its side (Figure 7A). The impact force was measured as the contact force between the impactor and the thorax. The shoulder deflection was defined as the change of distance between two nodes corresponding to the locations on humeri measured in the experiments (Figure 7B).

Table-Top Restraint-Like Loading Tests
Based on the investigation of the force-deflection response of the thorax under dynamic, non-impact, restraint-like loading at a non-injurious level, Kent et al. (2004) presented testing corridors subjected to single and double diagonal belts and hub and distributed loading on the anterior thorax. The limb amputated PMHS subjects were positioned supine on a table, and then compressed up to 20% of external chest depth for all four loading conditions. Four table-top FE models were set up to simulate the experimental loading cases (Figure 8). The pulleys were modeled using slip ring elements and 1-D belt elements to maintain the loading angles similar to experimental settings. The contact force between the back side of FBM and the table was output as the reaction force. The relative compression was characterized in terms of the ratio of the deflection of a point at the mid-sternum divided by initial distance between this point and the table along anteroposterior direction.

Model Performance Evaluation
To quantitatively evaluate the model performance, the outputs between model simulations and the experimental testing were compared by a widely used objective rating tool known as CORelation and Analysis or CORA (Gehre et al., 2009;Vavalle et al., 2013;Poulard et al., 2015). The CORA can evaluate the similarity between curves by the cross correlation method, which evaluates error according to phase shift, magnitude and curve shape to produce a relative score ranging from 0 (no correlation) to 1 (perfect match) (Miller et al., 2017;Decker et al., 2020). Supplementary Figures 3, 7, except the results of FSU flexion bending testing and rear hub impact for T-spine. The CORA ratings were collected in the Table 1.

Thoracic Intervertebral Joints
The simulation results of adjacent vertebrae to quasi-static external loads (force vs. disp., or moment vs. rotation) produced reasonable responses (see Supplementary Figure 3) compared to experimental data (Markolf, 1972). The average CORA score was 0.85, showing good correlation with the testing data. For dynamic flexion testing of upper-thoracic and mid-thoracic FSU,   ligament rupture was observed in the deformed configuration of the updated model, e.g., ligaments between T8 and T9 ( Figure 9A). The moment-angle curve generated by our updated model exhibited the toe-region followed by increasing stiffness that is characteristic of FSU joint response (Figures 9B,C). A region with sub-traumatic damage appeared until the peak moment before a significant drop in response, which matched the experimental results. The simulated peak moment before failure was 31.8 Nm for FSU T2-T4 and 48.3 N·m for FSU T7-T9, consistent with the experimental data. The average CORA score was 0.81, showing an overall reasonable correlation with the limited testing data. These results in the flexion simulation showed good agreement with soft tissue injuries and failure moment described in the observed experimental testing.
Following the FSU validation, a validation case involving a hub impact to the spine was simulated and a deformed configuration was shown in Figure 10A. The impact loading to the T-spine resulted in gross chest compression and forced the T-spine into extension. For all three impact speed tests, the impact force response produced by the model fell within the upper and lower bounds of the experimental data (Forman et al., 2015). For the 5.5 m/s impact, the peak chest deflection was close to the lower bound of the experimental data ( Figure 10B). For the spine extension angle, the simulation was close to the upper bound of the experiment data ( Figure 10C). Quantitatively, the CORA score for chest deflection was 0.94, and for the change in spine angle it was 0.89, showing a very good correlation with the testing data.

Costovertebral Joints
The results for all simulated costovertebral joints under each loading direction were summarized and provided in Supplementary Figure 4. Each of the cases had 2 to 3 momentangle curves from experimentation associated with it, and often there was considerable variation in the apparent stiffness.
Nevertheless, the CORA scores of most of the simulation results was 0.7 or greater, and the response was within the range of the experiments, except the ribs 2 under caudal motion. There were also some discrepancies for dorsal motion, which often produced stiffer responses likely due to contacts between the ribs and the transverse processes. It is unknown whether a similar effect was present in Duprey's tests (Duprey et al., 2010).

Clavicle and Sternal Modeling
The experiment testing reported reaction forces (along the direction of loading) against the deflection, and the peak bone surface strain against the forces. To evaluate the clavicle response, the simulation results were compared to the response obtained from three-point bending tests and axial compression tests, including the force versus the deflection and the peak cortical bone strain versus the force. It was noted that the response of the base model was softer when compared to the testing data. After a parametric study of the material properties, the elastic modulus of cortical clavicle bone was updated from 9 to 18 GPa, and the yield strength was increased from 0.08 to 0.16 GPa. The force-deflection response and peak strainforce response obtained from updated cortical bone model were plotted, which were reasonable compared to testing data (see Supplementary Figure 5).
For the three-point bending test of the sternum, the vertical reaction force was measured from the two supporting posts, which were associated with the proximal and distal ends of the sternum. The vertical force response obtained from the base model showed much higher peak force than the experiments. After performing parametric studies of material properties on cortical bone, the cortical sternum bone modulus was reduced from 14 to 4 GPa, and reducing yield strength from 0.09 to 0.035 GPa. The 2% MPS threshold was adopted to represent the strain based failure feature of cortical sternum bone. The vertical reaction force response on

Costal Cartilage Modeling
The experiment reported force-displacement results of the paired tests: perichondrium-intact and perichondrium-removed tests for each specimen. The peak forces (along anteroposterior direction) from each of the perichondrium-removed tests were normalized by dividing them by the peak forces of the matching perichondrium-intact tests. The simulation results of the model with optimized material parameters (the elastic modulus E = 12.5 MPa for the solid cartilage and E = 55 MPa for the fabric perichondrium shell layer) were compared to the results of the experimental data (see Supplementary Figure 7).
Compared to the FE model with perichondrium layer, the peak stress obtained from the model without perichondrium decreased about 45%, which was very similar to the 47% for the average ratio of measured peak force from the experiment.

Point Loading of the Eviscerated Ribcage
The predicted force-deflection responses from the ribcage model were compared with the experimental data (Figure 11). The rib 1 was overly stiff to a certain extent (i.e., the force was about 14% larger than the upper limit of the corridor when the maximum deflection was reached), and the rib 6 was close to the lower bound of the corridor. The CORA scores ranged from 0.81 to 0.91. Overall, the model responses agreed with the experimental corridors, which indicated that the ribcage stiffness was comparable to the tested three subjects.

Frontal Pendulum Impact
The sectional view of the thorax configuration under the frontal pendulum impact (25 ms) was provided in Figure 12A, which showed the ribcage and internal organs were compressed along the anteroposterior direction. The deflection and force responses predicted by the model were mostly within the range of the experimental corridors (Figures 12B,C). Quantitatively, the CORA score for chest deflection was 0.98, and for force response was 0.90, showing a very good correlation with the testing results.

Shoulder Pendulum Impact
In Figure 13, the model responses were compared with the experimental corridors, which were developed by Koh et al. (2005) based on biomechanical responses from four cadaveric studies. The overall deflection responses under different loading conditions were reasonable compared to the corridors, except for the 6.8 m/s unpadded impact, where the peak deflection predicted by the model was smaller to some extent compared to the range of the corridors. Compared to the corridor under each loading condition, the simulation predicted peak force appeared earlier and the value was lower than the corridor. The average CORA score of the chest deflection for the four loading cases was 0.90, and it was 0.85 for force response.

Table-Top Restraint-Like Loading Tests
In Figure 14, the simulated reaction force versus chest compression responses agreed well with the corresponding experimental corridors (Kent et al., 2004). The response curve produced by the double diagonal belts loading was slightly outside of the corridor bounds during some chest compression rates (e.g., less than 7% compression). However, the response curve produced by each of the other three loading conditions was within the corridor range. Quantitatively, the CORA ranged from 0.80 to 0.90, which indicated good correlation between simulation and testing results.

DISCUSSION
In this study, a hierarchical process was presented to evaluate, validate, and improve the musculoskeletal modeling of a midsized male thorax for impact scenarios. At the component level, the modeling enhancement and improvement for major bony structures of the thorax was introduced, including the thoracic intervertebral joints, costovertebral joints, clavicle, sternum and costal cartilages. At ribcage and full thorax levels, the thorax model with validated components was quantitatively evaluated by four types of experimental testing.

Body Parts Validation and Improvement
Most published thoracic spine models for impact analysis were part of thoracic or whole body models with oversimplified model features (Yang, 2018). For example, the GHBMC M50 v4.5 model (version released in 2016) used drum-like closed surfaces to represent intervertebral discs, and the spinous ligaments were not included in the model. In this work, the ligaments were modeled using available material property data from available literature, and the annulus fibrosus fiber lamina layers were created and embedded into the solid discs. Three levels of independent experimental tests have been conducted for validating the response of our enhanced T-spine model: adjacent vertebrae/segment, FSU, and full torso. For adjacent vertebrae in response to quasi-static loads (Markolf, 1972), the results of force vs. displacement and moment vs. rotation of T8-T9 showed good agreement with the experimental results.
However, it should be noted that the magnitudes of deformation or rotation under these quasi-static loads were very small (e.g., the maximum displacement in the reported shear experiments was less than 0.8 mm), which means slight variations between model setup and experiment were likely to cause discrepancies. The FSU flexion bending testing applied loading until soft tissue failure, which could describe the injury behavior of the upper and mid-thoracic spine at the segmental level. The initial toe region of the moment response and ligament rupture were well reproduced in the results (Figure 9) of our enhanced T-spine model with realistic IVD and non-linear ligaments. It was noted that the experimental study only included two specimens, which produced a wide range of response curves and failure moments. The blunt impacts on the back of thorax validated the chest deflection and spine extension angle of model, which showed the overall spine kinematic response under rear impact agreed well with the experimental results. For noninjurious impact testing cases (1.5 and 3.0 m/s), although the responses of deflection and change of spine extension angle were not presented here, they showed good agreement with the experimental results. The 5.5 m/s impact velocity was designed for potentially injurious situation. However, it lacked spinal injury in the experiment for the first two subjects with a nominal stroke distance of 85 mm. The impact stroke was increased to 150 mm for the final two subjects. Our simulation was conducted using 150 mm stroke distance. This explains why the simulated spine extension angle was at the upper limit of the data based on four subjects.  The simulation of the costovertebral joints under quasi-static ventral-dorsal flexion and cranial-caudal flexion predicted fair to good results compared to Duprey et al. (2010) test data, except the moment-angle produced by rib 2 caudal flexion. This could be attributed to the limited sample and often the samples themselves having wide range of response. Furthermore, the process of tuning the material properties of the ligaments associated with the costovertebral joint was complicated by the fact that many of the beams in the model, whether they represent the costotransverse, superior costotransverse, or lateral costotransverse ligaments, all couple and contribute to the moment-angle relationship for the 4 simulated directions. Stiffening up the beams to produce a better caudal response appears to also cause substantial deviations in the cranial, ventral, and dorsal relations. In a larger scope, multiple output responses of a computational model could be sensitive to the changes of a specific input parameter. it is important to ensure proper characterization of such parameters that have a significant impact on model outputs.
For clavicle and sternal validation and improvement, the material parameters were determined based on parametric studies of material properties. The adjusted material properties can compensate for missing model details or geometry deficiencies (Cronin, 2014), aiming to achieve better prediction of global kinematic response and component level injury. For clavicle under three-point bending and axial compression (Zhang et al., 2014), the predicted responses were within the range of experimental curves and the CORA scores were good, although the strain vs. force response was close to the upper bound of experimental data. Qualitatively, the trends of the simulated results of the three-point bending of sternum agreed well with the testing data (Kerrigan et al., 2010). However, the model only showed fair agreement with the experimental results according to the quantitative CORA scores, probably due to the large variability observed in testing data. To enhance the costal cartilage model, a layer of perichondrium with a fabric material was created to surround the solid mesh of the costal cartilage. Biomechanically, the perichondrium has been shown to contribute approximately half of the force response under cantilever-beam like bending , when comparing the normalized force-displacement curves produced by perichondrium-intact model and perichondrium-removed model. A similar distribution of reaction forces was achieved by tuning the Young's moduli for the costal cartilage and perichondrium. The results were within the testing range and reported experimental error.

Ribcage and Full Thorax Levels Validation
Once the presented body components were enhanced by adding more biofidelic details or adjusting material properties, and validated by adopted experimental loading cases, these components were integrated into the thorax in the FBM. To check the stability and biomechanical responses of thorax model, validations were conducted under a wide variety of loading scenarios at ribcage or full thorax level.
In this study, the validated sternum, costal cartilages, and thoracic spine were integrated into the thorax model. The geometry and material properties of the ribs in the current thorax model were not changed since they were already validated in previous studies Poulard et al., 2015;Guleyupoglu et al., 2017;Zaseck et al., 2018). The quasi-static point loading of the eviscerated ribcage  was utilized in this work to evaluate the ribcage model response with updated compoments. It should be noted that experimental corridors were developed based on only three denuded and eviscerated PMHS. The simulated force vs. deflection curves followed the expected trends established in the corridors, and were lower than the average experimental results, except in the case with point loading on the costochondral junction of rib level 1 (Rib1_CCJ). The frontal pendulum impact (Kroell et al., 1971) was simulated to evaluate the global chest response along the anteroposterior direction, including the chest deflection and impact force. The shoulder pendulum impact (Koh et al., 2005) was modeled to evaluate the overall response of the upper thorax along the lateral direction under four loading speed cases, including the deflection and impact force. Most of these pendulum impact cases were generally well correlated with the corresponding experimental validation corridors, except some local regions. This could be attributed to the modeling deficiencies of the shoulders (Xu et al., 2018), which lacked credible experimental testing data for validation thus far. Also, some details regarding how the corridors were produced were unknown, which could be a source of discrepancy between simulation and testing corridors. The predicted chest force vs. compression responses under the tabletop restraint-like loading tests (Kent et al., 2004) were within the upper and lower bounds of testing data. Quantitatively, the average CORA score was 0.88, which indicated the model responses closely agreed with the corridors.
In summary, to develop a credible thorax FE model, there is a need for a computational approach that can evaluate and validate model response for the sake of achieving high biofidelity. This work presented a hierarchical approach to validate and enhance thorax model biofidelity. Based on relevance and necessity to improve the current GHBMC thorax model, as well as the availability of experimental data, the aforementioned test cases were selected for validation at both component and whole thorax or torso levels. Overall, the model responses showed good or close agreement with the experimental data both qualitatively and quantitatively.

Limitations of This Study
Detailed human thorax body models can provide insights to investigate injury mechanisms and tolerances by qualitative assessment and quantitative evaluation of thoracic response at tissue level. However, similar to developing other human body models, there are some limitations of this work which may or may not affect the model biofidelity. First, there were a paucity of available experimental data or missing technical details regarding specific aspects of model verifications or validations. Unfortunately, most experimental data were not developed or presented for the purposes of computational modeling studies (Yang, 2018). Also, some experimental tests utilized for our model validation were conducted more than three or four decades ago. Some details related to the testing subjects or test conditions were missing or lacked adequate description. For instance, the frontal pendulum impact on chest (Kroell et al., 1971) were published in the 1970s.
Second, there were some model simplifications due to insufficient representation of anatomical details or material properties. For example, all ligaments related to the thoracic intervertebral joints and costovertebral joints were modeled using 1D elements, rather than a faithful representation of the 3D geometry. Also, the regional differences in cortical bone thickness were not considered in this thorax model, except for the ribs. Third, the model presented in the current study was not able to accurately predict local strains of the thoracic skeleton. At the component level validation, the strain-force response of the clavicle was compared between simulation and experimental results, since the data was available. However, reliable experimental data is needed to evaluate the strain distribution elsewhere within the whole ribcage. Even though a model can be validated against global impact responses, the local strains of most parts at tissue level are not able to be validated since cortical strain was relatively insensitive to the model inputs (Anderson et al., 2005) in a full body loading environment. While it was noted that the rib modeling and fractures of this presented thorax model has been validated in previous studies, a newer probabilistic method may be a better computational strategy in future investigations to predict rib fracture risk within FBM (Forman et al., 2012).
Finally, the current model was created based on a single midsized subject, which is not a real population-averaged model. The medical images for CAD data of this model were collected from a young (26-year-old) living male volunteer (Gayzik et al., 2011). Age-related factors affecting the thorax modeling such as morphologic and material characteristics were not investigated. The human ribcage morphology changes as it ages through adulthood. Also, due to the increased porosity and decreased mineralization associated with aging, the thinning of the cortical shell and changing bone material properties could change the biomechanics and injury tolerance of the thorax (Kent et al., 2005). For example, the tolerable sternal deflection level is much lower in the aging bony thorax. It is worth emphasizing that the modeling and validation work presented in this study were focused on biofidelity and capability of injury prediction, and further research is needed to validate the model for other applications (e.g., orthopedic biomechanics).

CONCLUSION
This work presented a comprehensive evaluation and validation procedure to hierarchically improve and enhance the MSK system of a mid-sized male thorax FE model. The material characteristics of relevant thoracic bony structures and ligaments were either improved by parametric study or determined based on available literature. To achieve high biofidelity, new model features were developed and incorporated into the chest components, such as thoracic intervertebral joints and costal cartilage. For the sake of accurately predicting injury risk at both component and whole thorax levels, a multi-level validation process was presented in this study to maximize the capability of FE model-predicted impact responses. The whole ribcage or thorax model was evaluated against a wide range of loads and different modes of loading. Compared with those corresponding experimental data sets, the CORA ratings of the model ranged from 0.76 to 0.99, indicating the model provides good predictions of the overall biomechanical responses of the whole thorax. This validated, highly biofidelic thorax model has been integrated into the newest generation of GHBMC models to represent the stateof-the-art computational human body models for crash injury prediction and prevention. Future studies will focus on tissuelevel injury predictions based on reliable experimental data, as well as extending the framework of accurate modeling and multi-level validation procedure for different population groups including age-and gender-based computational models.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
All authors contributed to the conception, design, and interpretation of results presented in this study, and are accountable for all aspects of this work.