Calibration of a Heterogeneous Brain Model Using a Subject-Specific Inverse Finite Element Approach

Central to the investigation of the biomechanics of traumatic brain injury (TBI) and the assessment of injury risk from head impact are finite element (FE) models of the human brain. However, many existing FE human brain models have been developed with simplified representations of the parenchyma, which may limit their applicability as an injury prediction tool. Recent advances in neuroimaging techniques and brain biomechanics provide new and necessary experimental data that can improve the biofidelity of FE brain models. In this study, the CAB-20MSym template model was developed, calibrated, and extensively verified. To implement material heterogeneity, a magnetic resonance elastography (MRE) template image was leveraged to define the relative stiffness gradient of the brain model. A multi-stage inverse FE (iFE) approach was used to calibrate the material parameters that defined the underlying non-linear deviatoric response by minimizing the error between model-predicted brain displacements and experimental displacement data. This process involved calibrating the infinitesimal shear modulus of the material using low-severity, low-deformation impact cases and the material non-linearity using high-severity, high-deformation cases from a dataset of in situ brain displacements obtained from cadaveric specimens. To minimize the geometric discrepancy between the FE models used in the iFE calibration and the cadaveric specimens from which the experimental data were obtained, subject-specific models of these cadaveric brain specimens were developed and used in the calibration process. Finally, the calibrated material parameters were extensively verified using independent brain displacement data from 33 rotational head impacts, spanning multiple loading directions (sagittal, coronal, axial), magnitudes (20–40 rad/s), durations (30–60 ms), and severity. Overall, the heterogeneous CAB-20MSym template model demonstrated good biofidelity with a mean overall CORA score of 0.63 ± 0.06 when compared to in situ brain displacement data. Strains predicted by the calibrated model under non-injurious rotational impacts in human volunteers (N = 6) also demonstrated similar biofidelity compared to in vivo measurements obtained from tagged magnetic resonance imaging studies. In addition to serving as an anatomically accurate model for further investigations of TBI biomechanics, the MRE-based framework for implementing material heterogeneity could serve as a foundation for incorporating subject-specific material properties in future models.

Introduction

Traumatic brain injury (TBI) is a significant source of injury, disability, and death. Recent epidemiological studies have estimated that TBI’s account for approximately one-third of all injury-related deaths in the United States (Faul and Coronado, 2015). In 2010, the Centers for Disease Control (CDC) estimated that TBIs resulted in 2.5 million emergency department (ED) visits (87%), hospitalizations (11%), and deaths (2%) (Centers for Disease Control and Prevention, 2015). Finite element (FE) models of the brain have rapidly become indispensable tools for investigating TBI mechanisms, assessing new protective technology, and developing injury risk criteria (Kleiven and von Holst, 2002; Giordano and Kleiven, 2014; Gabler et al., 2018). FE models of the brain are typically used to investigate the dynamic 3D deformation of the human brain under simulated head impacts relevant to sports, automotive crashes, and falls. While FE models have been instrumental to furthering our understanding of TBI biomechanics, many FE brain models have been developed to represent the brain as a simplified physical system, both in their representation of the anatomy and material properties, thus limiting their accuracy and utility in predicting deformations experienced due to head impacts.

Perhaps the most significant differences across brain models in the field relate to the constitutive laws and material parameters chosen to represent the material behavior of the simulated brain tissue (Jin et al., 2013; Dixit and Liu, 2017; Fahlstedt et al., 2021). At the simplest level, the brain is modeled as a single isotropic and homogeneous material (Kleiven and von Holst, 2002; Takhounts et al., 2008; Ji et al., 2015). However, unique material properties can be assigned to different parts of the brain, which typically represent tissue types with different cellular composition or segmented anatomical labels, and many FE models include differences in material properties between white and gray matter (Horga and Gilchrist, 2003; Kimpara et al., 2006; McAllister et al., 2012; Panzer et al., 2012; Mao et al., 2013; Miller et al., 2016, 2017). At the most complex level, the brain has been modeled as an anisotropic, heterogeneous structure by explicitly modeling axonal fiber tracts by embedding 1D elements in the brain mesh (Garimella et al., 2019; Hajiaghamemar et al., 2019; Wu et al., 2019). While the embedded axon approach provides a more biofidelic structural representation of the axonal tissue, which is known to exhibit anisotropy and regional variations in tissue material properties (Jin et al., 2013; Budday et al., 2015, 2017a, 2019; Weickenmeier et al., 2016), the embedding of 1D axonal tract elements can significantly increase computational cost. A potential alternative for modeling brain heterogeneity is to model the brain with material stiffness varying throughout the brain regardless of the tissue classification of each element, which can be obtained using magnetic resonance elastography (MRE).

Recently, MRE has been utilized to non-invasively measure in vivo material properties of the human brain in healthy volunteers (Weaver et al., 2012; Johnson et al., 2013a, b; Hiscox et al., 2016). In MRE, an external transducer (commonly a head pillow) mechanically vibrates the head (10–100 Hz) to induce micron-level displacements in the brain that can be measured and used to estimate elastic and viscous material properties throughout the brain with high spatial resolution (Hiscox et al., 2016). To date, MRE has been used to investigate global brain material properties with voxel-level resolution, regional variations in tissue stiffness, and material properties associated with brain pathology (Hiscox et al., 2016; Johnson and Telzer, 2018; Murphy et al., 2019). While stiffness measurements vary, most studies agree that the measured brain stiffness is dependent on the actuation frequency (due to viscoelasticity), and that white matter regions are stiffer than gray matter regions (Hiscox et al., 2016). However, since the measured properties are dependent on the actuation frequency and obtained from micron-level displacements of the brain, additional work is needed to apply the MRE-derived stiffness maps to a FE brain model.

The goal of this study is to calibrate and verify a heterogenous FE brain model by leveraging experimental datasets reporting (1) MRE-derived material properties (Hiscox et al., 2020), (2) high rate, in situ, brain displacements measured from human cadaveric specimens using sonomicrometry (Alshareef et al., 2020a), and (3) low rate, in vivo, brain strain measured from human volunteers using tagged magnetic resonance imaging (tMRI; Knutsen et al., 2020). Material parameters for the model were developed in three phases. First, stiffness data from an MRE template image (average of 134 subjects; Hiscox et al., 2020) was used to define the relative stiffness gradient throughout the brain model. Second, the linear stiffness parameter was calibrated using low-displacement cases from the Alshareef et al. (2020a) dataset and verified using brain strain data from the in vivo tMRI dataset. Finally, the non-linear stiffness parameter was calibrated using high-deformation cases from the Alshareef et al. (2020a) dataset. To verify that the calibrated material parameters were physically reasonable, a comprehensive verification was performed using the remaining rotational cases from these datasets that were not used for calibration. The response of the calibrated material used in the model was also compared to experimental in vitro material test data available in the literature (Jin et al., 2013) to ensure that the median stiffness response was within the range of experimental data. In addition to serving as a model for investigating TBI biomechanics, the technique for implementing MRE-derived heterogeneous material properties can be adapted to implement subject-specific material properties in future subject-specific brain models.

Materials and Methods

Development of CAB-20MSym Template Model

The mesh of the CAB-20MSym template model was chosen to represent the anatomy of the CAB-20MSym template image developed by Giudice et al. (2020). This template was constructed from T1-weighted MRI scans obtained from 20 young, healthy adult males (22 ± 3 years). Details regarding these images are provided elsewhere (Giudice et al., 2020; Reynier et al., 2020).

The CAB-20MSym template image was segmented to identify the brain parenchyma, peripheral cerebrospinal fluid (CSF), internal CSF, and ventricles. To generate the template model mesh, each 1 mm isotropic voxel in the segmentation image was directly converted into a cubic hexahedral (i.e., voxel) element and assigned to a part based on its segmentation label. This approach was selected as voxel meshes have the accuracy and stability benefits of hexahedral elements and can capture complex anatomical features at the native spatial resolution (in this case, 1 mm) of the MRI images used to construct them (Miller et al., 2016; Ghajari et al., 2017). To include the sagittal sinus, falx cerebri, and tentorium cerebelli, the voxels surrounding these regions were manually delineated and two-dimensional shell elements were generated at the corresponding mid-surfaces. Finally, a layer of rigid shell elements surrounding the outermost surface of the peripheral CSF part was generated to represent the dura, which was assumed to be rigidly connected to the inner surface of the skull (Miller et al., 2016). The selection of the numerical implementation approach was informed by an analytical review of the numerical methods utilized by the TBI modeling community performed by Giudice et al. (2019b) and is summarized in the Supplementary Material. All interfaces were continuous and connected through shared nodes.

The material properties for the CSF (peripheral and internal), ventricles, skull, sagittal sinus, falx, and tentorium were adapted from previous brain models (Takhounts et al., 2008; Mao et al., 2013; Miller et al., 2016). To account for the nominal stiffness provided by the trabeculae and bridging vessels located within the subarachnoid space, the peripheral CSF was modeled using a linear viscoelastic material with very low stiffness (G₀ = 0.5 kPa; G_∞ = 0.1 kPa). As these properties do not exist for CSF, an elastic fluid (bulk modulus, K = 2.1 GPa) was assigned to the ventricle and internal CSF parts. The sagittal sinus, falx, and tentorium were modeled as elastic materials (Young’s modulus, E = 31.5 MPa; Poisson’s ratio, ν = 0.45). Finally, the skull was modeled as rigid to allow the implementation of 6 degree-of-freedom head kinematic boundary conditions for all analyses (Gabler et al., 2016). Further details regarding these material properties are available in the Supplementary Material. The material implementation of the brain parenchyma is described in Sections “Implementation of Brain Heterogeneity” and “Constitutive Modeling of Brain Parenchyma.”

Implementation of Brain Heterogeneity

The implementation of brain material heterogeneity in the CAB-20MSym template model was derived from the MRE134 template image (Hiscox et al., 2020). This template was constructed using MRE data from 134 healthy, young adults (18–35 years, 78F/56M) using common MRE acquisition and data processing protocols (Hiscox et al., 2020). To adapt the MRE134 template, originally defined in MNI152 space with 2 mm isotropic voxels, it was first non-linearly registered to the CAB-20MSym template space using ANTs non-linear registration (Avants et al., 2008). To eliminate stiffness measurements potentially influenced by numerical artifacts or edge effects, stiffness values beyond the 98th percentile were excluded. In the CAB-20MSym template model, CSF spaces were modeled using CSF-specific constitutive models, however, in the MRE134 template image, CSF spaces were not differentiated from the brain tissue. Therefore, stiffness values below the 15th percentile, which corresponded to the approximate stiffness in the MRE134 template in CSF areas (approximately 0–1.5 kPa; Hiscox et al., 2020), were excluded. From this truncated distribution, the stiffness value of each voxel was normalized by the median stiffness and binned into 10 groups. Voxels that had original stiffnesses below the 15th and above the 98th percentile were assigned to the lowest and largest normalized stiffness bins, respectively. This process yielded a normalized stiffness label image where voxels were categorized by their relative stiffness, and not according to an anatomical segmentation label. In doing so, this approach accounts for stiffness variations present within tissue groups, such as white and gray matter. Normal variations in the material properties of various tissue groups are reported in the literature (Hiscox et al., 2020).

Constitutive Modeling of Brain Parenchyma

The 10 parts comprising the brain parenchyma, using the binned groups from the previous section, were modeled using a quasi-linear viscoelastic (QLV) model (Fung, 1993). In Fung’s QLV theory, it is assumed that the response of a material can be separated into a normalized function of time only, g(t), and an elastic function of strain only, T^e(ε) (Fung, 1993). For a QLV model, the stress relaxation function, R(ε,t), is:

$$R(\mathrm{\epsilon },t)=g(t)\cdot T^e(\mathrm{\epsilon })$$

In this study, the instantaneous elastic response function, T^e(ε), was derived from an Ogden strain energy density (Ogden and Hill, 1972).

$$W({\mathrm{\lambda }}_1,{\mathrm{\lambda }}_2,{\mathrm{\lambda }}_3)=\frac{\mathrm{\mu }}{\mathrm{\alpha }}\left({\mathrm{\lambda }}_1^{\mathrm{\alpha }}+{\mathrm{\lambda }}_2^{\mathrm{\alpha }}+{\mathrm{\lambda }}_3^{\mathrm{\alpha }}-3\right)$$

Where λ_j are the three principal stretches, μ is the shear modulus, and α is a unitless non-linearity coefficient. The infinitesimal shear modulus (i.e., initial slope of the non-linear shear stress-strain curve), μ₀, can be obtained as a function of the material parameters.

$${\mathrm{\mu }}_0=\frac{1}{2}\mathrm{\mu }\mathrm{\alpha }$$

Finally, four Prony terms (N = 4) were included in the reduced relaxation function.

$$g\left(t\right)=g_{\mathrm{\infty }}+\sum _{i=1}^N{g}_i{e}^{-{\mathrm{\beta }}_{i}t}$$

$${\mathrm{\beta }}_i=\frac{1}{{\mathrm{\tau }}_i}$$

$$g_{\mathrm{\infty }}+\sum _{i=1}^N{g}_i=1$$

Where, g_∞ and g_i are normalized coefficients associated with the long-term response and each time constant, τ_i. The density, ρ, and Poisson’s ratio, ν, of the brain parenchyma material were 1.123 × 10^–6 kg/mm³ (Miller et al., 2016) and 0.499999, respectively.

The reduced relaxation function was fit to experimental tan(δ) data obtained from studies that characterized the viscoelastic properties of brain tissue over wide ranges of input frequencies (Fallenstein et al., 1969; Shuck and Advani, 1972; Arbogast et al., 1997; Bilston et al., 1997, 2001; Arbogast and Margulies, 1998; Brands, 2000; Darvish and Crandall, 2001; Lippert et al., 2004; Nicolle et al., 2004; Hrapko et al., 2006; Shen et al., 2006; Garo et al., 2007). The tan(δ) response represents material damping as a function of frequency, and is independent of the material stiffness, making it suitable for calibrating the reduced relaxation function. Optimization of g_i and τ_i was performed using a least squares optimization. The Ogden parameters defining the median deviatoric response of the brain, to which the brain parenchyma heterogeneity was applied relative to, was calibrated using an iFE calibration scheme and is described in the following sections.

Calibration Objective and Approach

To reduce the number of optimized parameters and simplify the calibration process, it was assumed that the damping of the brain and hyperelastic non-linearity were homogeneous throughout the brain. As such, the reduced relaxation function parameters were assigned to all brain parenchyma parts. Therefore, the objective of this inverse FE approach was to calibrate two parameters:

1. The median stiffness of the brain, μ_med, that defined the material heterogeneity.

2. The non-linearity coefficient, α, that defined the hyperelastic non-linearity.

Material parameters were calibrated using an iFE approach in which parameters were optimized to minimize the error between model and experimental results. In this case, a subset of the impact cases in the in situ brain displacement dataset (Alshareef et al., 2018, 2020a, b) were simulated using subject-specific models and the error between nodal displacements and the corresponding experimental displacements were minimized. In this dataset, in situ brain deformation was measured at discrete locations using sonomicrometry sensors (i.e., “receivers”), distributed throughout the brain parenchyma of cadaveric head-neck specimens. Three subject-specific models were generated to represent the anatomies of subjects SONO-896, SONO-900, and SONO-904 (Table 1; Alshareef et al., 2020a). Further details regarding this experimental dataset and the corresponding model setup are provided in the Supplementary Materials.

TABLE 1

Table 1. Specimen information for subjects used to assess deformation response (Alshareef et al., 2020a).

To minimize the likelihood of obtaining a non-unique solution, μ_med and α were optimized independently by leveraging specific cases in the in situ brain displacement database. To optimize μ_med a low severity case (Axial 20 rad/s, 60 ms) was used to optimize the median infinitesimal quasi-static shear stiffness (a function of μ_med). This case was selected as preliminary simulations indicated that the deformations induced in this case were not sensitive to the material non-linearity coefficient. A higher severity case (Coronal 40 rad/s, 30 ms) was used to optimize α while maintaining the optimized infinitesimal quasi-static shear stiffness. This coronal rotation case was selected as preliminary simulations indicated that the deformations predicted in this case were sensitive to the material non-linearity coefficient and avoided calibrating the material parameters under a single loading direction. Finally, the calibrated heterogeneous material parameters were verified using the remaining cases in the in situ brain deformation database, as well as low severity in vivo brain strain data obtained using the in vivo tMRI database (Knutsen et al., 2020). This final verification included 39 simulations; 11 simulations for each SONO subject-specific model (not including case used for calibration) and 6 simulations for the tMRI cases. Details regarding these experimental datasets and the respective model setups are provided in the Supplementary Materials. As a final check, the material response was assessed in tension, compression, and simple shear loading over strain rates of 0.5, 5, and 30 1/s and compared to experimental data (Jin et al., 2013).

The material calibration process is summarized below:

1. Optimize μ_med under low severity rotational impacts.

2. Intermediate verification to assess calibrated μ_med under independent low severity cases.

3. Optimize α under higher severity rotational impacts.

4. Comprehensive verification of calibrated material parameters

All models used in this study were subject-specific, constructed to represent the anatomies of the subjects from which the experimental data were collected. As such, the deformation response of each subject-specific model was compared to the experimental data obtained from that subject only. This was crucial to eliminate geometric effects from the calibration and verification processes, ensuring that discrepancies between simulation and experimental results, for each subject, were due to variations in material parameters only. All subject-specific models were generated using registration-based morphing (RBM), which is a non-linear morphing technique that was developed specifically for generating subject-specific models of the brain by leveraging image registration transformations (Giudice et al., 2020). Metrics of morphing accuracy and element quality for all subject-specific models generated in this study are provided in the Supplementary Materials.

Objective Rating and Optimization Strategy

CORrelation and Analysis (CORA) scores were used to quantify the error between the nodal displacements in the model and the receiver displacements in the experimental dataset. Since experimental displacements were 3D (i.e., x, y, and z displacement time-histories for each receiver), a composite score (cCORA) was calculated to obtain a single objective rating for each receiver. cCORA was computed as the weighted average of the CORA scores in each orthogonal direction, weighted by the relative magnitude of the experimental signal in each direction (Giudice et al., 2019a). Finally, the overall score for each rotational case was computed as the weighted average of all receiver cCORA scores, weighted by the experimental maximum resultant displacement for each receiver (wcCORA). In this study, the default CORA parameters were used (Gehre et al., 2009).

$$wcCORA=\sum _i^N{\mathrm{\alpha }}_i\times cCOR{A}_i$$

$${\mathrm{\alpha }}_i=\frac{{\mathrm{\beta }}_{\text{i}}}{{\sum }_i^N{\mathrm{\beta }}_{\text{i}}}$$

Where, β_i is the experimental maximum resultant displacement for the i^th receiver and N is the number of receivers for each rotation case. For the calibration of the Ogden material parameters, wcCORA was used to quantify the error between the nodal displacements and the experimental receiver displacements for each rotation simulation. To obtain a set of parameters that best represented the overall response of the three subjects used for calibration, a joint optimization was performed where the goal was to maximize the mean wcCORA score for subjects SONO-896, SONO-900, and SONO-904.

$$f(x)=\frac{wcCOR{A}_{896}+wcCOR{A}_{900}+wcCOR{A}_{904}}{3}$$

A golden ratio search algorithm was used to identify the maximum mean overall wcCORA as a function of the material parameter being optimized. In the first iteration, a series of simulations were run to identify the bounds for the golden search algorithm. For example, when calibrating α, three simulations with α = 2, 6, and 10 were run to determine the bounds for calibration. The parameters investigated in each subsequent iteration, i (x_i,1 and x_i,2) were determined by the golden ratio search algorithm.

$$x_{i,1}=a_i+\left(1-\mathrm{\phi }\right)\left(b_i-a\right)$$

$$x_{i,2}=a_i+\mathrm{\phi }\left(b_i-a\right)$$

$$\mathrm{\phi }=\frac{\sqrt{5}-1}{2}$$

Where a_i and b_i are the lower and upper bounds for iteration, i. The bounds were updated based on the mean wcCORA score for each parameter investigated [f(x₁) and f(x₂)] in the previous iteration (i – 1).

$$\mathrm{If}\mathrm{f}\left({\mathrm{x}}_1\right)>\mathrm{f}\left({\mathrm{x}}_2\right)\mathrm{then}{\mathrm{a}}_{\mathrm{i}}={\mathrm{a}}_{(\mathrm{i}-1)}\mathrm{and}{\mathrm{b}}_{\mathrm{i}}={\mathrm{x}}_{(\mathrm{i}-1),2}$$

$$\mathrm{If}\mathrm{f}\left({\mathrm{x}}_1\right)<\mathrm{f}\left({\mathrm{x}}_2\right)\mathrm{then}{\mathrm{a}}_{\mathrm{i}}={\mathrm{x}}_{(\mathrm{i}-1),1}\mathrm{and}{\mathrm{b}}_{\mathrm{i}}={\mathrm{b}}_{(\mathrm{i}-1)}$$

This process was repeated until the termination criteria was satisfied. The parameter that had the greatest mean wcCORA was selected as the calibrated value. The first iteration and termination criteria for the calibration of μ_med and α are shown in Table 2. These values were selected based on parameters reported in the literature for human brain tissue (Miller and Chinzei, 2002; Nicolle et al., 2004; Franceschini et al., 2006; Kleiven, 2007; Kaster et al., 2011; Moran et al., 2014; Budday et al., 2017a).

TABLE 2

Table 2. First iteration parameters and termination criteria for calibrated material parameters.

Calibration and Verification of Median Shear Modulus

In the first step of the material calibration process, μ_med was optimized by simulating the Axial 20 rad/s, 60 ms case (Z: 20–60). All simulations were run for 200 ms. Pilot simulations indicated that the predicted deformations in this loading case were not sensitive to material non-linearity that is governed by α. Therefore, in this optimization the shear response of the material was constrained to the linear response of a Neo-Hookean solid (α = 2), where the shear modulus defined in the Ogden model is identical to the infinitesimal shear modulus of the material (μ₀).

Since the model incorporates material heterogeneity relative to the median stiffness of the brain, the stiffness of each of the 10 parts was defined relative to the median stiffness defined for any given optimization iteration.

$${\mathrm{\mu }}_i={\mathrm{\gamma }}_i{\mathrm{\mu }}_{med}={\mathrm{\gamma }}_i{\mathrm{\mu }}_{0,med}$$

Where γ_i is the relative stiffness of the ith part and varies between 0.53 and 1.53. All stiffness values reported in this study refer to quasi-static parameters (i.e., μ_med and α represent the quasi-static response of the material).

To verify the calibrated median shear stiffness (μ_0,med = μ_med) the remaining 20 rad/s, 60 ms cases in the coronal (X: 20–60) and sagittal (Y: 20–60) directions were simulated and compared to each specimen’s experimental data using wcCORA. Verification of the strain response was performed by simulating in vivo brain strain experiments (Knutsen et al., 2020). In this dataset, three subjects (tM-3978, tM-4838, and tM-6176) were subjected to sagittal rotations (ω_max = 1.4–1.6 rad/s) and three subjects (tM-3978, tM-7126, and tM-9475) were subjected to axial rotations (ω_max = 4–5.4 rad/s) of the head. This verification step was performed to ensure that the calibrated μ_0,med was physically meaningful and resulted in biofidelic deformation predictions under independent low-deformation test cases. This also provided reassurance that the final set of calibrated material parameters were unique, given the model’s sensitivity to both μ and α.

To compare the predicted strain response of each subject-specific model, the maximum principal strain (MPS) of each element was computed and mapped to the corresponding voxel in the subject image. As a global metric of strain, the 95th percentile MPS value (MPS-95) predicted by the model was compared to the equivalent experimental measures. In addition, the volume fraction of elements exceeding 2% strain was compared between the model and experimental data. These volume fractions were computed globally, as well as regionally for the cerebral gray and white matter and the cerebellum. These regions were identified using a segmentation image provided in the tMRI database. MPS-95 was also computed for the elements/voxels located in these regions.

Calibration and Verification of Non-linear Coefficient

In the previous steps, the median shear stiffness, which was equivalent to the median infinitesimal shear modulus (μ_0,med = μ_med) since the shear response was assumed to be linear (α = 2), was determined and verified. In this step, the non-linear coefficient (α) was calibrated by simulating the Coronal 40 rad/s, 30 ms rotation case. All simulations were run for 200 ms. To preserve the previously calibrated μ_0,med, μ_med was also adjusted such that μ_0,med was equivalent to the value determined in the first step of the optimization procedure. Therefore, for the i^th part in the heterogeneous CAB-20MSym template model the shear modulus in the Ogden constitutive model, μ_i, was defined as a function of μ_med, α, and γ_i.

$${\mathrm{\mu }}_i=\frac{2{\mathrm{\gamma }}_i{\mathrm{\mu }}_{0,med}}{\mathrm{\alpha }}$$

To verify the calibrated heterogeneous material, the remaining 11 rotation cases for each subject were simulated (36 total simulations including case used for calibration), and the nodal displacements of the calibrated subject-specific models were compared to the corresponding experimental brain displacement data using wcCORA. To assess the calibrated model performance relative to other state-of-the-art and widely used FE brain models, wcCORA values obtained in this study were compared to those reported for the Global Human Body Models Consortium (GHBMC) brain model (Mao et al., 2013) and the UVA embedded axon model (UVA-EAM) (Wu et al., 2019, 2020, 2021). These models were morphed to the anatomy of the three subjects using surface-based morphing (Wu et al., 2019) and simulated under identical boundary conditions, resulting in 36 simulations per model.

The set of tMRI simulations used to verify the calibrated median infinitesimal shear stiffness were also simulated using the final calibrated heterogeneous model to ensure that the incorporation of material non-linearity did not influence the model predictions. Since the brain strains in these simulations were low (less than 6% MPS), material non-linearity was not expected to influence the model results.

Finally, the calibrated Ogden material response was compared to experimental material test data. To verify the optimized parameters, the complex modulus of the derived model was compared to the complex moduli reported in the rheological characterization dataset. Furthermore, a series of single element (1 mm × 1 mm × 1 mm) simulations were run with the fit material parameters to verify the response of the constitutive model as well as the Ogden QLV implementation in LS-Dyna. These single element simulations were run in tension, compression, and simple shear at loading rates of 0.5, 5, and 30 1/s to 50% engineering strain and the results were compared to average response corridors constructed from the material characterization data in the literature (Jin et al., 2013).

Results

Development of CAB-20MSym Template Model

The CAB-20MSym template model had an intracranial volume of 1439 cm³ and approximately 1.6 million elements, 1.5 million nodes, and 16 parts (Figure 1). Note that Figure 1 depicts the external surfaces of each of these parts and is not representative of the entire volume of each part. All interfaces (e.g., falx-brain) were continuous and defined using shared nodes.

FIGURE 1

Figure 1. Depiction of the parts of the CAB-20MSym mesh. Brain parenchmya parts are subdivided based on their stiffness relative to the median stiffness (μ/μ_med) of the MRE134 template. The volume of each solid part is indicated in parentheses. The skull is shown in gray.

Implementation of Brain Heterogeneity

The original and truncated stiffness distributions from the MRE134 template, mapped to CAB-20MSym space, are shown in Figure 2. In the truncated distribution, the median stiffness was 2.53 kPa and the mean ± standard deviation was 2.37 ± 0.99 kPa. To convert the MRE134 template image into the CAB-20MSym template model, the truncated distribution was grouped into 10 bins, which were used to classify the normalized stiffness of each voxel in the CAB-20MSym model. The 10 relative stiffness groups were 0.53, 0.64, 0.75, 0.86, 0.97, 1.08, 1.20, 1.31, 1.42, and 1.53, which corresponded to the bin centers (Figure 2). The MRE134 template, normalized stiffness image, and normalized stiffness label image are shown in Figure 3.

FIGURE 2

Figure 2. (Left) original (blue) and truncated (orange) MRE134 stiffness distributions. Stiffness values in the original distribution were truncated between the 15th and 98th percentile stiffness values. (Right) grouped normalized stiffness values (N = 10 bins).

FIGURE 3

Figure 3. Magnetic resonance elastography template (MRE134), normalized stiffness image, and normalized stiffness label image (used to construct CAB-20MSym template model). Each of the 10 bins in the normalized stiffness label corresponded to an individual part in the template model.

Constitutive Modeling of the Brain Parenchyma

The reduced relaxation function parameters were obtained by fitting g_i and β_i (i = 4) to match experimental tan(δ) data (Figure 4). The fit reduced relaxation parameters are shown in Table 3. These viscoelastic parameters were applied homogeneously throughout the entire brain model.

FIGURE 4

Figure 4. tan(δ) (left) and complex modulus (right) of the calibrated material parameters, demonstrating the median complex modulus of the brain (purple) as well as the complex modulus for the softest (green) and stiffest (blue) relative stiffness parts. Note that tan(δ) is a function of damping only and not tissue stiffness. Experimental data were obtained from the literature (see Section “Constitutive Modeling of Brain Parenchyma” for list of references).

TABLE 3

Table 3. Prony series parameters for the brain parenchyma.

Calibration and Verification of Median Shear Modulus

The median shear stiffness of the material was calibrated using the axial 20 rad/s, 60 ms rotation case in the in situ brain displacement database. The calibrated value of μ_med was 1.125 kPa. Since this initial constitutive model constrained the non-linear coefficient to 2 (representing a Neo-Hookean solid), μ_med was equivalent to the infinitesimal quasi-static shear modulus, μ_0,med. Six iterations (12 simulations) were required to satisfy the termination criteria specified in Table 2, resulting in an optimal mean wcCORA of 0.62 (Figure 5). The calibrated value was within the range of other Ogden rubber infinitesimal quasi-static shear moduli reported in the literature for the brain (Miller and Chinzei, 2002; Nicolle et al., 2004; Franceschini et al., 2006; Kleiven, 2007; Kaster et al., 2011; Moran et al., 2014; Budday et al., 2017a), which varied from 0.27 to 1.49 (mean ± standard deviation = 0.92 ± 0.38).

FIGURE 5

Figure 5. Calibration results for μ_med (left) and α (right). Parameters were calibrated to maximize the mean wcCORA score (black) across subjects SONO-896 (blue), SONO-900 (red), and SONO-904 (yellow). wcCORA for each subject was obtained by comparing nodal displacements from corresponding subject-specific models and experimental data.

Using the jointly calibrated μ_med = μ_0,med value, the Coronal and Sagittal 20 rad/s, 60 ms cases for the same specimens were simulated to verify the subject-specific model responses. The wcCORA scores for the coronal, sagittal, and axial simulations were between 0.66–0.72, 0.50–0.67, and 0.58–0.69, respectively (Table 4). In general, wcCORA scores were greatest in the coronal rotation, and across the three subjects, the SONO-900 model yielded the highest wcCORA scores using the jointly calibrated μ_med value.

TABLE 4

Table 4. wcCORA scores for μ_0,med joint optimization and verification.

A series of tMRI experiments were simulated using subject-specific models to verify the calibrated value of μ_med = μ_0,med. The strain response of subject tM-9475 during the axial rotation is shown in Figure 6. Qualitatively, the tM-9475 subject-specific model demonstrated a similar deformation pattern compared to the experimental data. In general, MPS was largest in the cortex and smallest in the cerebellum (Figure 6). Additionally, in both simulation and experiment, larger strains were observed at the apexes of the ventricles (i.e., frontal and occipital horns) which propagated through the cerebrum. This can be seen in slices z = 20–40 in Figure 6. However, in the experimental data cortical strain was consistently asymmetric, with the largest strains observed in the right hemisphere of the brain, whereas in the simulation the MPS distribution was relatively symmetric. These trends were similar for all subjects tested in the axial rotation protocol.

FIGURE 6

Figure 6. Comparison of maximum principal strain distribution for subject tM-9475 under axial rotation (top) and tM-4838 under sagittal rotation (bottom).

These trends were salient in the MPS-95 and strain volume fraction metrics used to quantify brain deformation globally and regionally. For subject tM-9475, global and regional MPS-95 values were similar to the experimental data, with absolute differences between 0.0025 and 0.0084 strain. However, while regional MPS-95 values were similar between the model and experiment, larger differences in the volume fraction of voxels that exceeded 2% strain was observed, particularly in the cerebellum (simulation: 0.15; experiment: 0.33) and gray matter (simulation: 0.34; experiment: 0.61). Figures showing MPS-95 and volume fraction results for the axial and sagittal tMRI cases are included in the Supplementary Materials.

The strain response of subject tM-4838 during the sagittal rotation is shown in Figure 6. Qualitatively, the tM-4838 subject-specific model demonstrated a similar MPS pattern compared to the experimental data with the largest MPS values observed at the periphery of the brain, including the cortex, base of the cerebrum, and anterior surface of the brainstem, and lowest in the midbrain (Figure 6). However, strains in the cerebellum and brainstem were larger in the experimental data than predicted by the simulation. This is particularly evident in slices x = −15–5 in Figure 6. In general, global and regional MPS-95 metrics were similar between the subject-specific model and experimental data with absolute differences in MPS-95 between 0.0003 and 0.0067 for subject tM-4838. The volume fraction of voxels that exceeded 2% strain was also similar, except for the cerebellum, where the experimental volume fraction was 0.16 compared to only 0.06 in the simulation. These trends were similar for all subjects tested in the sagittal rotation protocol.

Calibration and Verification of Non-linear Coefficient

The non-linearity parameter of the Ogden constitutive model, α, was calibrated using the Coronal 40 rad/s, 30 ms rotation case in the in situ brain deformation database. The calibrated value of α was 6.67, which yielded a mean wcCORA of 0.67 across subjects SONO-896, SONO-900, and SONO-904 (Figure 5). Eight iterations (12 simulations) were required to satisfy the termination criteria specified in Table 2. The final calibrated material parameters for each brain parenchyma part (grouped by stiffness relative to the median stiffness) are shown in Table 5. The viscoelastic parameters for the brain parts are shown in Table 3.

TABLE 5

Table 5. Calibrated Ogden material parameters for each brain parenchyma part, grouped by stiffness relative to the median value.

To verify the calibrated material parameters, all 12 impact cases for subjects SONO-896, SONO-900, and SONO-904 were simulated and assessed using wcCORA. The three subject-specific models demonstrated good biofidelity with a mean wcCORA of 0.63 ± 0.06 for all 36 simulations (range: 0.50–0.74). Across the three subjects, SONO-900 demonstrated the highest wcCORA scores (0.68 ± 0.03) and wcCORA was similar for SONO-896 (0.61 ± 0.07) and SONO-904 (0.61 ± 0.06). In general, wcCORA scores were highest in the coronal impacts (0.68 ± 0.04) and similar in the sagittal (0.61 ± 0.06) and axial (0.61 ± 0.06) impacts. A summary of the wcCORA scores for all three subjects is shown in Figure 7. Exemplary nodal displacement-time histories for subject SONO-904 in the 40 rad/s, 60 ms cases are shown in the Supplementary Materials.

FIGURE 7

Figure 7. Comparison of wcCORA scores for the calibrated heterogeneous model developed in this study, the GHBMC brain model (Mao et al., 2013), and the anisotropic UVA-EAM model (Wu et al., 2019).

The tMRI simulations were run with the final calibrated material model to ensure that the addition of material non-linearity did not influence the strain prediction. As expected, the brain strain response in these simulations was dominated by the infinitesimal shear stiffness, which was optimized in the first stage of material calibration. The incorporation of material non-linearity had negligible effect on the MPS distribution since this second optimization was performed while maintaining the infinitesimal shear modulus previously calibrated.

The complex modulus, as a function of frequency, for the median stiffness response and softest (μ/μ_med = 0.53) and stiffest (μ/μ_med = 1.53) relative stiffness parts are shown in Figure 4. All three responses were within the range of experimental data. Finally, the median, softest (μ/μ_med = 0.53), and stiffest (μ/μ_med = 1.53) material responses were assessed in tension, compression, and simple shear loading at strain rates of 0.5, 5, and 30 1/s and compared to experimental data (Figure 8; Jin et al., 2013). At the highest strain rate (30 1/s) and intermediate (5 1/s), the three material responses were like the experimental corridors for all three loading modes, simultaneously. At the lowest strain rate (0.5 1/s), the three material responses were similar to the experimental in compression but were stiffer in tension and shear.

FIGURE 8

Figure 8. Calibrated material response in tension, compression, and simple shear at rates of 0.5, 5, and 30 1/s compared to experimental material characterization data (Jin et al., 2013). The median (purple), softest (green), and stiffest (blue) response of the brain are shown.

Discussion

In this study the heterogeneous CAB-20MSym template model was developed, calibrated, and extensively evaluated using a comprehensive set of experimental data that included measurements of material properties using in vivo MRE experiments (Hiscox et al., 2020), in situ brain displacement measured using sonomicrometry (Alshareef et al., 2020a), in vivo brain strain measured using tagged MRI (Knutsen et al., 2020), and in vitro material response (Jin et al., 2013; Meaney et al., 2014). In all assessments of model biofidelity, the CAB-20MSym template model demonstrated a high fidelity to the experimental data which represented a spectrum of TBI severity, ranging from non-injurious (tMRI cases) to moderate-to-severe TBI (sonomicrometry cases).

Our approach to implement heterogeneity using MRE data was selected for several reasons. This technique implemented material heterogeneity without significantly increasing the computational cost. Other studies have implemented heterogeneity using embedded beam elements to explicitly model the structural contributions of axonal fiber tracts (Garimella and Kraft, 2017; Garimella et al., 2019; Wu et al., 2019). However, this technique has been reported to increase the computational cost by a factor of 2.4. For the models developed to calibrate the material parameters, this would have increased the computational time for each simulation from approximately 15 to 36 h (simulating 200 ms of response). As such, the computational time required to run the battery of the 36 in situ brain deformation simulations (Alshareef et al., 2020a) would have increased from approximately 540 to 1300 h.

The relative stiffness approach was also chosen with consideration of effects specific to the MRE data. The MRE stiffness measurements we used were determined at a single frequency of 50 Hz, though brain tissue properties from MRE have been shown to be dependent on the actuation frequency due to the viscoelasticity of the brain (Sack et al., 2013). Thus, the specific shear stiffness values recovered at this frequency may not be most relevant for a TBI model. Furthermore, MRE measures tissue stiffness under micron-level displacements. Since the brain is a highly non-linear material, these stiffness results may not be applicable at the finite levels of deformation associated with TBI. Using relative stiffness minimizes any contributions these factors may have had on the absolute stiffness measurements.

There are also several limitations associated with the MRE-based approach. Firstly, although brain heterogeneity was represented as a function of tissue stiffness and not tissue type, the implementation was still isotropic as brain MRE measurements typically assume material isotropy, though methods for extracting anisotropic properties are being developed (Romano et al., 2012; Tweten et al., 2015; McGarry et al., 2020; Smith et al., 2020). Secondly, it was assumed that the relative stiffnesses obtained from the micron-level displacements used in MRE were linear with strain, and that the heterogeneity observed in MRE was applicable for higher levels of deformation. While the validity of these assumptions is unclear, the verification of the brain response with these material parameters over a wide range of rotational loading conditions minimizes their significance.

An inverse FE approach was used to calibrate the material parameters of the heterogeneous CAB-20MSym template model using subject-specific models. However, three challenges faced by any iFE calibration problem are overfitting, obtaining a unique solution, and computational cost. The likelihood of overfitting and obtaining non-unique solutions in an optimization problem is typically increased by using an excessively complicated model (i.e., optimizing many material parameters) and an imbalance between the amount of data used to fit the model and the data used to validate the model (i.e., excessive training data). Furthermore, since one simulation is required for each test case, the computational cost of calibrating a material model increases significantly with the number of parameters included. Given that each simulation in this study took approximately 15 h to run, a robust and efficient approach to material calibration was required.

In the heterogeneous CAB-20MSym template model, each brain parenchyma material requires two parameters to describe bulk properties (density, ρ and Poisson’s ratio, ν), two parameters to describe deviatoric elastic properties (μ and α), and eight parameters to describe viscoelastic properties (g₁–g₄ and β₁–β₄), for a total of 12 parameters per material. Without further reduction, a full material calibration would require optimizing 120 parameters, which would likely result in non-unique solutions and the computational cost would be prohibitive. In this study, several assumptions were made to reduce the number of material parameters included in the optimization problem. First, all brain tissues were assumed to have the same density and Poisson’s ratio, which is consistent with the TBI biomechanics literature and is common practice in FE brain modeling (Takhounts et al., 2008; Mao et al., 2013; Wu et al., 2019; Alshareef et al., 2020a). This reduced the number of material parameters from 120 to 100. By assuming homogeneous damping and non-linearity, the total number of variables calibrated in this study was reduced to 11 (μ for each material, and a single non-linearity coefficient,α). For the stiffness variables, each material was defined relative to the median value using experimental MRE data, further reducing the total number of calibrated parameters to 2 (μ_med and α). While non-linearity and viscoelastic properties of the brain have been shown to vary spatially throughout the brain (Johnson et al., 2016; Budday et al., 2017a, b, 2019; Hiscox et al., 2018), these are more challenging measurements and the extent to which viscoelasticity and elastic non-linearity vary remains an open question. These disparities are likely attributed to differences in experimental protocols, including tissue harvest sites, tissue hydration, and loading conditions (Budday et al., 2019). The membranes implemented in the model (falx and tentorium) and CSF parts were not included in the optimization, as a preliminary sensitivity study indicated that the influences of the stiffness and damping of these structures on the deformation response was negligible, compared to the brain parenchyma parts.

We used a combined optimization to calibrate the median shear stiffness and non-linear coefficient, as opposed to an individual calibration for each of the three subjects. The objective of the study was to obtain a singular set of optimal material parameters to be used in the CAB-20MSym template model. Subject-specific material properties could have been combined (e.g., averaged), however, due to the highly non-linear nature of this optimization problem this may not have yielded an optimal solution for the CAB-20MSym template model. Nonetheless, the calibration process indicated that the three subjects likely had different underlying material properties (Figure 5). For example, the estimated optimal median shear moduli for subjects SONO-896, SONO-900, and SONO-904 were approximately 0.70, 1.18, and 2.05 kPa, respectively. A simple average of these values would have resulted in an “optimal” shear modulus of 1.31 kPa, compared to 1.13 kPa determined by the optimization which maximized wcCORA. An important limitation was that only three elderly (57–67 years), cadaveric subjects were used in the optimization process due to the availability of experimental brain deformation data. However, the calibrated model demonstrated a biofidelic response when verified using the in vivo tMRI dataset, which included younger, living subjects between 21 and 42 years of age (Knutsen et al., 2020). Nonetheless, the calibrated material parameters may not be representative of the general population and this calibration should be repeated once brain deformation data for more specimens is available. It should also be noted that since these material parameters were obtained using an inverse FE process, they are model specific and cannot be arbitrarily applied to other models (Giudice et al., 2019b). However, given that the calibrated materials have been verified using experimental material characterization data, it is very likely that they represent the underlying material response, and not the contributions of the numerical implementation. Nonetheless, if implemented in another model, the deformation response should be thoroughly investigated prior to model deployment.

The final calibrated material model was verified by its response to experimental tissue data from multiple sources. The complex modulus of the calibrated materials and the response in tension, compression, and simple shear were compared to experimental material characterization data for cortical tissue. In both assessments, the median, softest, and stiffest material responses in the calibrated CAB-20MSym template model were similar the experimental data, verifying that the material parameters obtained from the iFE calibration were physically meaningful and representative of the underlying tissue. In these comparisons, the full range of material response was compared to the experimental data to account for the fact that the experimental data may not necessarily represent the median response of the brain.

The biofidelity of the calibrated CAB-20MSym template model was assessed using brain deformation data from the experiments conducted by Alshareef et al. (2018, 2020a) and Knutsen et al. (2020). Collectively, these datasets encapsulated various magnitudes (2.5–40 rad/s), durations (30–60 ms), and directions of rotational loading (coronal, sagittal, axial). The model demonstrated biofidelic responses under the test conditions used by Alshareef et al. Interestingly, differences in wcCORA scores across subjects reflected the subject-specific relationships between material properties and wcCORA scores in the optimization process. For example, wcCORA scores for SONO-900 were generally greater than those for SONO-896. In the material calibration process, the estimated optimal material properties for SONO-900 were similar to the joint optimal value (μ_0,med = 1.18 kPa), whereas the estimated optimal material properties for SONO-896 were softer (μ_0,med ≈ 0.7 kPa). The opposite was true for SONO-904 (μ_0,med ≈ 2.0 kPa). These differences in calibrated parameters likely reflect variability in material properties between these specimens. While variation in material parameters is to be expected, future work is required to quantify this variability across larger populations. This may be an important factor as subject-specific models of the brain become more prominent in research and clinical applications. Moving forward, the MRE-based framework used to implement material heterogeneity in this study will be adapted to incorporate subject-specific material properties in future models that consider the entire subject-specific brain structure. While the CAB-20MSym template model demonstrated a biofidelic response under the loading conditions investigated in this study, additional verification is required if the model is exercised under loading conditions that deviate from those used in the current study (e.g., blast loading or micron-level harmonic displacements).

While a high-resolution measurement of brain strain at injurious loading conditions is not available, assessing the biofidelity of a brain model’s strain prediction using low severity, non-injurious loading conditions can improve the confidence of a model’s strain response. In general, the CAB-20MSym demonstrated good biofidelity under the axial and sagittal tMRI cases compared to the experimental data. However, there were some discrepancies between the model and experimental results, particularly in the volume fraction exceeding 2% strain. These differences could have been attributed to several factors. Firstly, in the experimental tMRI data set, only resultant head cradle kinematics were recorded. Therefore, to simulate these experiments it was assumed that the head was perfectly coupled to the head cradle and that the applied rotational head kinematics were perfectly uniaxial (either axial or sagittal). While the head was tightly coupled to the head cradle, it is possible that there were slight discrepancies between the head kinematics in the experiments and simulations. Secondly, the brain-skull interface was implemented by modeling the CSF layer between the brain and skull. While many approaches have been investigated to model this interface (Wang et al., 2018), the relative motion between the brain and skull is not well characterized and there is no consensus on best modeling practices for this interface. Therefore, discrepancies between the relative skull-brain motion in the simulations and experiments could have contributed to the observed differences in volume fraction and deformation fields, especially in these low deformation impacts. Furthermore, the brainstem of the CAB-20MSym template model was truncated at the foramen magnum due to a lack of MRE data in the inferior portions of the brainstem and superior portions of the spinal cord. This may have attributed to some of the observed differences, especially in the inferior regions of the brain. However, since the deformations induced in these experiments were small, it is possible that these factors had an exaggerated effect on the predicted strains, and it is not clear how these effects translate to larger deformation cases. Finally, while tMRI is a well-established imaging technique, it can be susceptible to experimental error (approximately 0.7% strain, introduced during filtering or interpolation), which could affect these low strain measurements and the resulting 2% strain volume fractions (Gomez et al., 2019). Nonetheless, the overall biofidelity of the CAB-20MSym template model’s strain prediction under these low severity loading conditions was satisfactory.

Summary

In this study, the CAB-20MSym template model was developed, calibrated, and extensively verified over a wide range of rotational head kinematic loading conditions. This model utilized a computationally efficient approach for incorporating material heterogeneity that leveraged data from a MRE template image that represented the average brain stiffness of 134 healthy adult subjects. Overall, the developed model demonstrated a biofidelic response for both nodal displacement and element strain metrics. Moving forward this template model will serve as the foundation of the registration-based morphing pipeline developed by Giudice et al. (2020) and can serve as an anatomically accurate model for further investigations of TBI mechanisms and to aid the development of novel protective equipment and safety countermeasures. Furthermore, the framework for implementing material heterogeneity using MRE data can be adapted to incorporate subject-specific material properties in future models of the brain.

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found here: NKI-RS: http://fcon_1000.projects.nitrc.org/indi/enhanced/; MRE-134 Template: https://github.com/mechneurolab/mre134; In Situ Brain Deformation: https://www.nhtsa.gov/research-data/databases-and-software.

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.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors gratefully acknowledge the Achievement Rewards for Collegiate Scientists (ARCS) Foundation and the UVA Brain Institute for supporting this research.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fbioe.2021.664268/full#supplementary-material

References

Alshareef, A., Giudice, J. S., Forman, J., Salzar, R. S., and Panzer, M. B. (2018). A novel method for quantifying human in situ whole brain deformation under rotational loading using sonomicrometry. J. Neurotrauma 35, 780–789. doi: 10.1089/neu.2017.5362

PubMed Abstract | CrossRef Full Text | Google Scholar

Alshareef, A., Giudice, J. S., Forman, J., Shedd, D. F., Reynier, K. A., Wu, T., et al. (2020a). Biomechanics of the human brain during dynamic rotation of the head. J. Neurotrauma 37, 1546–1555. doi: 10.1089/neu.2019.6847

PubMed Abstract | CrossRef Full Text | Google Scholar

Alshareef, A., Giudice, J. S., Forman, J., Shedd, D. F., Wu, T., Reynier, K. A., et al. (2020b). Application of trilateration and Kalman filtering algorithms to track dynamic brain deformation using sonomicrometry. Biomed. Signal Process. Control 56:10. doi: 10.1016/j.bspc.2019.101691

CrossRef Full Text | Google Scholar

Arbogast, K. B., and Margulies, S. S. (1998). Material characterization of the brainstem from oscillatory shear tests. J. Biomech. 31, 801–807. doi: 10.1016/S0021-9290(98)00068-2

CrossRef Full Text | Google Scholar

Arbogast, K. B., Thibault, K. L., Pinheiro, B. S., Winey, K. I., and Margulies, S. S. (1997). A high-frequency shear device for testing soft biological tissues. J. Biomech. 30, 757–759. doi: 10.1016/s0021-9290(97)00023-7

CrossRef Full Text | Google Scholar

Avants, B., Epstein, C., Grossman, M., and Gee, J. (2008). Symmetric diffeomorphic image registration with cross-correlation: evaluating automated labeling of elderly and neurodegenerative brain. Med. Image Anal. 12, 26–41. doi: 10.1016/j.media.2007.06.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Bilston, L. E., Liu, Z., and Phan-Thien, N. (1997). Linear viscoelastic properties of bovine brain tissue in shear. Biorheology 34, 377–385. doi: 10.1016/s0006-355x(98)00022-5

CrossRef Full Text | Google Scholar

Bilston, L. E., Liu, Z., and Phan-Thien, N. (2001). Large strain behaviour of brain tissue in shear: some experimental data and differential constitutive model. Biorheology 38, 335–345.

Google Scholar

Brands, D. A. (2000). The large shear strain dynamic behavior of in-vitro porcine brain tissue and a silicone gel model material. Stapp Car Crash J. 44, 249–260.

Google Scholar

Budday, S., Nay, R., de Rooij, R., Steinmann, P., Wyrobek, T., Ovaert, T. C., et al. (2015). Mechanical properties of gray and white matter brain tissue by indentation. J. Mech. Behav. Biomed. Mater. 46, 318–330. doi: 10.1016/j.jmbbm.2015.02.024

PubMed Abstract | CrossRef Full Text | Google Scholar

Budday, S., Ovaert, T. C., Holzapfel, G. A., Steinmann, P., and Kuhl, E. (2019). Fifty shades of brain: a review on the mechanical testing and modeling of brain tissue. Arch. Comput. Methods Eng. 27, 1187–1230. doi: 10.1007/s11831-019-09352-w

CrossRef Full Text | Google Scholar

Budday, S., Sommer, G., Birkl, C., Langkammer, C., Haybaeck, J., Kohnert, J., et al. (2017a). Mechanical characterization of human brain tissue. Acta Biomater. 48, 319–340. doi: 10.1016/j.actbio.2016.10.036

PubMed Abstract | CrossRef Full Text | Google Scholar

Budday, S., Sommer, G., Haybaeck, J., Steinmann, P., Holzapfel, G. A., and Kuhl, E. (2017b). Rheological characterization of human brain tissue. Acta Biomater. 60, 315–329. doi: 10.1016/j.actbio.2017.06.024

PubMed Abstract | CrossRef Full Text | Google Scholar

Centers for Disease Control and Prevention (2015). “Report to congress on traumatic brain injury in the united states: epidemiology and rehabilitation,” in National Center for Injury Prevention and Control. Presented at the Division of Unintentional Injury Prevention, (Atlanta, GA: Centers for Disease Control and Prevention), 1–72. Available online at: https://www.cdc.gov/traumaticbraininjury/pdf/TBI_Report_to_Congress_Epi_and_Rehab-a.pdf

Google Scholar

Darvish, K. K., and Crandall, J. R. (2001). Nonlinear viscoelastic effects in oscillatory shear deformation of brain tissue. Med. Eng. Phys. 23, 633–645. doi: 10.1016/S1350-4533(01)00101-1

CrossRef Full Text | Google Scholar

Dixit, P., and Liu, G. R. (2017). A review on recent development of finite element models for head injury simulations. Arch. Comput. Methods Eng. 24, 979–1031. doi: 10.1007/s11831-016-9196-x

CrossRef Full Text | Google Scholar

Fahlstedt, M., Abayazid, F., Panzer, M. B., Trotta, A., Zhao, W., Ghajari, M., et al. (2021). Ranking and rating bicycle helmet safety performance in oblique impacts using eight different brain injury models. Ann. Biomed. Eng. 49, 1097–1109. doi: 10.1007/s10439-020-02703-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Fallenstein, G. T., Hulce, V. D., and Melvin, J. W. (1969). Dynamic mechanical properties of human brain tissue. J. Biomechan. 2, 217–226. doi: 10.1016/0021-9290(69)90079-7

CrossRef Full Text | Google Scholar

Faul, M., and Coronado, V. (2015). Epidemiology of Traumatic Brain Injury, in: Handbook of Clinical Neurology. Amsterdam: Elsevier, 3–13.

Google Scholar

Franceschini, G., Bigoni, D., Regitnig, P., and Holzapfel, G. A. (2006). Brain tissue deforms similarly to filled elastomers and follows consolidation theory. J. Mech. Phys. Solids 54, 2592–2620. doi: 10.1016/j.jmps.2006.05.004

CrossRef Full Text | Google Scholar

Fung, Y. (1993). Biomechanics: Mechanical Properties of Living Tissues. Berlin: Springer Science & Business Media.

Google Scholar

Gabler, L. F., Crandall, J. R., and Panzer, M. B. (2016). Assessment of kinematic brain injury metrics for predicting strain responses in diverse automotive impact conditions. Ann. Biomed. Eng. 44, 3705–3718. doi: 10.1007/s10439-016-1697-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Gabler, L. F., Crandall, J. R., and Panzer, M. B. (2018). Development of a second-order system for rapid estimation of maximum brain strain. Ann. Biomed. Eng. 47, 1971–1981. doi: 10.1007/s10439-018-02179-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Garimella, H. T., and Kraft, R. H. (2017). Modeling the mechanics of axonal fiber tracts using the embedded finite element method: axonal fiber mechanics using the embedded element method. Int. J. Numer. Methods Biomed. Eng. 33:e2823. doi: 10.1002/cnm.2823

PubMed Abstract | CrossRef Full Text | Google Scholar

Garimella, H. T., Menghani, R. R., Gerber, J. I., Sridhar, S., and Kraft, R. H. (2019). Embedded finite elements for modeling axonal injury. Ann. Biomed. Eng. 47, 1889–1907. doi: 10.1007/s10439-018-02166-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Garo, A., Hrapko, M., Van Dommelen, J. A. W., and Peters, G. W. M. (2007). Towards a reliable characterisation of the mechanical behaviour of brain tissue: the effects of post-mortem time and sample preparation. Biorheology 44, 51–58.

Google Scholar

Gehre, C., Gades, H., and Wernicke, P. (2009). “Objective rating of signals using test and simulation responses,” in Proceedings of the 21st International Technical Conference on the Enhanced Safety of Vehicles, (Stuttgart).

Google Scholar

Ghajari, M., Hellyer, P. J., and Sharp, D. J. (2017). Computational modelling of traumatic brain injury predicts the location of chronic traumatic encephalopathy pathology. Brain 140, 333–343. doi: 10.1093/brain/aww317

PubMed Abstract | CrossRef Full Text | Google Scholar

Giordano, C., and Kleiven, S. (2014). Evaluation of axonal strain as a predictor for mild traumatic brain injuries using finite element modeling. Stapp. Car Crash J. 58:29.

Google Scholar

Giudice, J. S., Alshareef, A., Wu, T., Gancayco, C. A., Reynier, K. A., Tustison, N. J., et al. (2020). An image registration-based morphing technique for generating subject-specific brain finite element models. Ann. Biomed. Eng. 48, 2412–2424. doi: 10.1007/s10439-020-02584-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Giudice, J. S., Park, G., Kong, K., Bailey, A., Kent, R., and Panzer, M. B. (2019a). Development of open-source dummy and impactor models for the assessment of american football helmet finite element models. Ann. Biomed. Eng. 47, 464–474. doi: 10.1007/s10439-018-02155-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Giudice, J. S., Zeng, W., Wu, T., Alshareef, A., Shedd, D. F., and Panzer, M. B. (2019b). An analytical review of the numerical methods used for finite element modeling of traumatic brain injury. Ann. Biomed. Eng. 47, 1855–1872. doi: 10.1007/s10439-018-02161-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Gomez, A. D., Knutsen, A. K., Xing, F., Lu, Y.-C., Chan, D., Pham, D. L., et al. (2019). 3-D measurements of acceleration-induced brain deformation via harmonic phase analysis and finite-element models. IEEE Trans. Biomed. Eng. 66, 1456–1467. doi: 10.1109/TBME.2018.2874591

PubMed Abstract | CrossRef Full Text | Google Scholar

Hajiaghamemar, M., Wu, T., Panzer, M. B., and Margulies, S. S. (2019). Embedded axonal fiber tracts improve finite element model predictions of traumatic brain injury. Biomech. Model. Mechanobiol. 19, 1109–1130. doi: 10.1007/s10237-019-01273-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Hiscox, L. V., Johnson, C. L., McGarry, M. D. J., Perrins, M., Littlejohn, A., van Beek, E. J. R., et al. (2018). High-resolution magnetic resonance elastography reveals differences in subcortical gray matter viscoelasticity between young and healthy older adults. Neurobiol. Aging 65, 158–167. doi: 10.1016/j.neurobiolaging.2018.01.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Hiscox, L. V., McGarry, M. D., Schwarb, H., Van Houten, E. E., Pohlig, R. T., Roberts, N., et al. (2020). Standard-space atlas of the viscoelastic properties of the human brain. Hum. Brain Mapp. 41, 5282–5300. doi: 10.1002/hbm.25192

PubMed Abstract | CrossRef Full Text | Google Scholar

Hiscox, L. V., Johnson, C. L., Barnhill, E., McGarry, M. D., Huston, J. III, Van Beek, E. J., et al. (2016) MRE of the human brain: technique, findings and clinical applications. Phys. Med. Biol. 61:R401.

Google Scholar

Horgan, T. J., and Gilchrist, M. D. (2003). The creation of three-dimensional finite element models for simulating head impact biomechanics. Int. J. Crashworthiness 8, 353–366. doi: 10.1533/ijcr.2003.0243

CrossRef Full Text | Google Scholar

Hrapko, M., Van Dommelen, J. A. W., Peters, G. W. M., and Wismans, J. (2006). The mechanical behaviour of brain tissue: large strain response and constitutive modelling. Biorheology 43, 623–636.

Google Scholar

Ji, S., Zhao, W., Ford, J. C., Beckwith, J. G., Bolander, R. P., Greenwald, R. M., et al. (2015). Group-Wise evaluation and comparison of white matter fiber strain and maximum principal strain in sports-related concussion. J. Neurotrauma 32, 441–454. doi: 10.1089/neu.2013.3268

PubMed Abstract | CrossRef Full Text | Google Scholar

Jin, X., Zhu, F., Mao, H., Shen, M., and Yang, K. H. (2013). A comprehensive experimental study on material properties of human brain tissue. J. Biomech. 46, 2795–2801. doi: 10.1016/j.jbiomech.2013.09.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, C. L., McGarry, M. D., Gharibans, A. A., Weaver, J. B., Paulsen, K. D., Wang, H., et al. (2013a). Local mechanical properties of white matter structures in the human brain. Neuroimage 79, 145–152. doi: 10.1016/j.neuroimage.2013.04.089

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, C. L., McGarry, M. D., Van Houten, E. E., Weaver, J. B., Paulsen, K. D., Sutton, B. P., et al. (2013b). Magnetic resonance elastography of the brain using multishot spiral readouts with self-navigated motion correction. Magn. Reson. Med. 70, 404–412. doi: 10.1002/mrm.24473

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, C. L., Schwarb, H., McGarry, M. D. J., Anderson, A. T., Huesmann, G. R., Sutton, B. P., et al. (2016). Viscoelasticity of subcortical gray matter structures. Hum. Brain Mapp. 37, 4221–4233. doi: 10.1002/hbm.23314

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, C. L., and Telzer, E. H. (2018). Magnetic resonance elastography for examining developmental changes in the mechanical properties of the brain. Dev. Cogn. Neurosci. 33, 176–181. doi: 10.1016/j.dcn.2017.08.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Kaster, T., Sack, I., and Samani, A. (2011). Measurement of the hyperelastic properties of ex vivo brain tissue slices. J. Biomech. 44, 1158–1163. doi: 10.1016/j.jbiomech.2011.01.019

PubMed Abstract | CrossRef Full Text | Google Scholar

Kimpara, H., Nakahira, Y., Iwamoto, M., Miki, K., Ichihara, K., and Kawano, S. (2006). Investigation of anteroposterior head-neck responses during severe frontal impacts using a brain-spinal cord complex FE model. Stapp. Car Crash J. 50:509.

Google Scholar

Kleiven, S. (2007). Predictors for traumatic brain injuries evaluated through accident reconstructions. Stapp. Car Crash J. 51, 81–114.

Google Scholar

Kleiven, S., and von Holst, H. (2002). Consequences of head size following trauma to the human head. J. Biomech. 35, 153–160. doi: 10.1016/s0021-9290(01)00202-0

CrossRef Full Text | Google Scholar

Knutsen, A. K., Gomez, A. D., Gangolli, M., Wang, W.-T., Chan, D., Lu, Y.-C., et al. (2020). In vivo estimates of axonal stretch and 3D brain deformation during mild head impact. Brain Multiphysics 1:100015. doi: 10.1016/j.brain.2020.100015

PubMed Abstract | CrossRef Full Text | Google Scholar

Lippert, S. A., Rang, E. M., and Grimm, M. J. (2004). The high frequency properties of brain tissue. Biorheology 41, 681–691.

Google Scholar

Mao, H., Zhang, L., Jiang, B., Genthikatti, V. V., Jin, X., Zhu, F., et al. (2013). Development of a finite element human head model partially validated with thirty five experimental cases. J. Biomech. Eng. 135:111002. doi: 10.1115/1.4025101

CrossRef Full Text | Google Scholar

McAllister, T. W., Ford, J. C., Ji, S., Beckwith, J. G., Flashman, L. A., Paulsen, K., et al. (2012). Maximum principal strain and strain rate associated with concussion diagnosis correlates with changes in corpus callosum white matter indices. Ann. Biomed. Eng. 40, 127–140. doi: 10.1007/s10439-011-0402-6

PubMed Abstract | CrossRef Full Text | Google Scholar

McGarry, M. D., Van Houten, E., Guertler, C., Okamoto, R. J., Smith, D. R., Sowinski, D. R., et al. (2020). A heterogenous, time harmonic, nearly incompressible transverse isotropic finite element brain simulation platform for MR elastography. Phys. Med. Biol. doi: 10.1088/1361-6560/ab9a84 [Epub ahead of print].

CrossRef Full Text | PubMed Abstract | Google Scholar

Meaney, D. F., Morrison, B., and Bass, C. D. (2014). The mechanics of traumatic brain injury: a review of what we know and what we need to know for reducing its societal burden. J. Biomech. Eng. 136:021008.

Google Scholar

Miller, K., and Chinzei, K. (2002). Mechanical properties of brain tissue in tension. J. Biomech. 35, 483–490. doi: 10.1016/S0021-9290(01)00234-2

CrossRef Full Text | Google Scholar

Miller, L. E., Urban, J. E., and Stitzel, J. D. (2016). Development and validation of an atlas-based finite element brain model. Biomech. Model. Mechanobiol. 15, 1201–1214. doi: 10.1007/s10237-015-0754-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Miller, L. E., Urban, J. E., and Stitzel, J. D. (2017). Validation performance comparison for finite element models of the human brain. Comput. Methods Biomech. Biomed. Eng. 20, 1273–1288. doi: 10.1080/10255842.2017.1340462

PubMed Abstract | CrossRef Full Text | Google Scholar

Moran, R., Smith, J. H., and García, J. J. (2014). Fitted hyperelastic parameters for human brain tissue from reported tension, compression, and shear tests. J. Biomech. 47, 3762–3766. doi: 10.1016/j.jbiomech.2014.09.030

PubMed Abstract | CrossRef Full Text | Google Scholar

Murphy, M. C., Huston, J. III, and Ehman, R. L. (2019). MR elastography of the brain and its application in neurological diseases. NeuroImage 187, 176–183. doi: 10.1016/j.neuroimage.2017.10.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Nicolle, S., Lounis, M., and Willinger, R. (2004). Shear properties of brain tissue over a frequency range relevant for automotive impact situations: new experimental results. Stapp. Car Crash J. 48, 239–258.

Google Scholar

Ogden, R. W., and Hill, R. (1972). Large deformation isotropic elasticity – on the correlation of theory and experiment for incompressible rubberlike solids. Proc. R. Soc. Lond. Math. Phys. Sci. 326, 565–584. doi: 10.1098/rspa.1972.0026

CrossRef Full Text | Google Scholar

Panzer, M. B., Myers, B. S., Capehart, B. P., and Bass, C. R. (2012). Development of a finite element model for blast brain injury and the effects of CSF cavitation. Ann. Biomed. Eng. 40, 1530–1544. doi: 10.1007/s10439-012-0519-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Reynier, K. A., Alshareef, A., Sanchez, E. J., Shedd, D. F., Walton, S. R., Erdman, N. K., et al. (2020). The effect of muscle activation on head kinematics during non-injurious head impacts in human subjects. Ann. Biomed. Eng. 48, 2751–2762. doi: 10.1007/s10439-020-02609-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Romano, A., Scheel, M., Hirsch, S., Braun, J., and Sack, I. (2012). In vivo waveguide elastography of white matter tracts in the human brain. Magn. Reson. Med. 68, 1410–1422. doi: 10.1002/mrm.24141

PubMed Abstract | CrossRef Full Text | Google Scholar

Sack, I., Jöhrens, K., Würfel, J., and Braun, J. (2013). Structure-sensitive elastography: on the viscoelastic powerlaw behavior of in vivo human tissue in health and disease. Soft. Matter. 9, 5672–5680. doi: 10.1039/c3sm50552a

CrossRef Full Text | Google Scholar

Shen, F., Tay, T. E., Li, J. Z., Nigen, S., Lee, P. V. S., and Chan, H. K. (2006). Modified Bilston nonlinear viscoelastic model for finite element head injury studies. J. Biomech Eng. 128, 797–801. doi: 10.1115/1.2264393

CrossRef Full Text | Google Scholar

Shuck, L. Z., and Advani, S. H. (1972). Rheological response of human brain tissue in shear. J. Basic Eng. 94, 905–911. doi: 10.1115/1.3425588

CrossRef Full Text | Google Scholar

Smith, D. R., Guertler, C. A., Okamoto, R. J., Romano, A. J., Bayly, P. V., and Johnson, C. L. (2020). Multi-Excitation magnetic resonance elastography of the brain: wave propagation in anisotropic white matter. J. Biomech. Eng. 142:071005.

Google Scholar

Takhounts, E. G., Ridella, S. A., Hasija, V., Tannous, R. E., Campbell, J. Q., Malone, D., et al. (2008). Investigation of traumatic brain injuries using the next generation of simulated injury monitor (SIMon) finite element head model. Stapp. Car Crash J. 52, 1–31. doi: 10.1155/2015/837585

PubMed Abstract | CrossRef Full Text | Google Scholar

Tweten, D. J., Okamoto, R. J., Schmidt, J. L., Garbow, J. R., and Bayly, P. V. (2015). Estimation of material parameters from slow and fast shear waves in an incompressible, transversely isotropic material. J. Biomech. 48, 4002–4009. doi: 10.1016/j.jbiomech.2015.09.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, F., Han, Y., Wang, B., Peng, Q., Huang, X., Miller, K., et al. (2018). Prediction of brain deformations and risk of traumatic brain injury due to closed-head impact: quantitative analysis of the effects of boundary conditions and brain tissue constitutive model. Biomechan. Model. Mechanobiol. 17, 1165–1185. doi: 10.1007/s10237-018-1021-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Weaver, J. B., Pattison, A. J., McGarry, M. D., Perreard, I. M., Swienckowski, J. G., Eskey, C. J., et al. (2012). Brain mechanical property measurement using MRE with intrinsic activation. Phys. Med. Biol. 57:7275. doi: 10.1088/0031-9155/57/22/7275

CrossRef Full Text | Google Scholar

Weickenmeier, J., de Rooij, R., Budday, S., Steinmann, P., Ovaert, T. C., and Kuhl, E. (2016). Brain stiffness increases with myelin content. Acta Biomater. 42, 265–272. doi: 10.1016/j.actbio.2016.07.040

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, T., Alshareef, A., Giudice, J. S., and Panzer, M. B. (2019). Explicit modeling of white matter axonal fiber tracts in a finite element brain model. Ann. Biomed. Eng. 47, 1908–1922. doi: 10.1007/s10439-019-02239-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, T., Antona-Makoshi, J., Alshareef, A., Giudice, J. S., and Panzer, M. B. (2020). Investigation of cross-species scaling methods for traumatic brain injury using finite element analysis. J. Neurotrauma 37, 410–422. doi: 10.1089/neu.2019.6576

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, T., Hajiaghamemar, M., Giudice, J. S., Alshareef, A., Margulies, S., and Panzer, M. B. (2021). Evaluation of tissue-level brain injury metrics using species-specific simulations. J. Neurotrauma Neu. 2020:7445. doi: 10.1089/neu.2020.7445

PubMed Abstract | CrossRef Full Text | Google Scholar