Modeling and simulation of brain herniation caused by subdural hematoma using finite-element method

. Brain shift and herniation are important signs of increased intracranial pressure (ICP) caused by hematomas or other types of intracranial mass. We propose a novel finite-element model that can be deformed in response to increased ICP. The half sphere model of the brain is partially divided into two compartments by the intact mid-sagittal plane, allowing subfalcine herniation. A 40 mm circle in the center of its equatorial plane allows transtentorial herniation. We perform a single load step, structural static analysis, simulating a left-sided subdural hematoma (SDH) compressing the cerebral hemispheres from the outer surface of the left hemisphere. Subfalcine and transtentorial brain herniations are reproduced and visualized. The Poisson’s ratio represents the tightness of the brain and the pressure load repre sents the ICP. There is a linear relationship between maximal deformation and the pressure load. The maximal deformation at the basal circumference and that at the basal midline closely resembles the maximal thickness of the SDH and the midline shift. We have developed a simple finite-element model that can simulate brain shift and herniation caused by pressure loads exerted on its surface by a mass. The experimental results correlate well with clinical observation on patients with acute and chronic SDH.


Introduction
Traumatic brain injury (TBI) is an important cause of mortality and severe disability worldwide. Pathophysiologically, TBI can be classified under the categories of primary and secondary injury, which are considered as stages in the evolution of head trauma [1]. By definition, "primary" head injury occurs at impact and may involve neural or vascular elements of the brain. Although surgeons can remove some intracranial hematoma, there is no treatment that can reverse primary injury to brain parenchyma and prevention is the only thing that helps. On the other hand, "secondary" injury is the result of physiological processes that are initiated following impact or soon thereafter. Caused by a Corresponding author : d95548001@ntu.edu.tw systemic or intracranial conditions, secondary injury may aggravate primary brain injury, but it can be ameliorated by taking specific measures.
After trauma, the edematous and deformed brain can produce increased intracranial pressure (ICP), resulting in brain shift, which is an important manifestation of secondary injury. The injured cerebral hemisphere may herniate through the subfalcine space (SFS) to the other side, resulting in subfalcine herniation (SFH), appearing as midline shift (MLS) on brain images such as computed tomography (CT) studies [2]. With further elevation of the ICP, the medial part of the cerebrum may herniate through the tentorial incisura (TI), resulting in transtentorial herniation (TTH) and brain stem compression. On images, one can find changes in the shape of the brain stem or effacement of the cisternal spaces around it [2].
Biomechanics is devoted to the analysis, measurement and modeling of the effects which are measured quantitatively under various mechanical loading situations in biological systems. The special discipline of biomechanics which is concerned with injury caused by mechanical interaction is denoted as trauma biomechanics [3]. We have proposed methods that can automatically measure the magnitude of MLS with symmetry-based and landmark-based methods [4][5][6]. To our knowledge, there has been no biomechanical model dedicated to simulating brain herniation in the literature.
Being the most important tool in trauma biomechanics, the finite-element (FE) method has been widely used in determining head injury criteria for motor vehicle accident and sports-related injuries. In addition to cadaver experiments based on global head kinematics, complex FE models of the head composed of different tissue types or complex geometry have been proposed to better describe regional response of the brain to primary injury [7][8][9].
The goal of this study is to develop the simplest 3-dimensional FE model that can be applied to quantitatively simulate brain deformation, specifically brain herniation, caused by an intracranial mass with increased ICP. The resulting brain deformation predicted by our model is then evaluated and compared to our previous results derived from real patients.

Finite-element model specification
We have constructed a simple model of the supratentorial brain [10]. The model assumes the shape of the upper half of a sphere with a diameter of 160 mm and a volume of 1072 ml. The outer surface of this model simulates the convexity surfaces of the cerebral hemispheres, while the irregular basal surface of the cerebrum is simplified and represented by a flat equatorial plane in this model. A 40 mm circle is removed from the center of the basal surface to simulate the midbrain, the upper part of the brain stem occupying the TI.
To model the partially connected cerebral hemispheres, our model is divided into two compartments by the intact mid-sagittal plane (iMSP). Each cerebral "hemisphere", which actually assumes the shape of quarter sphere, has three surfaces: the convexity, the basal and the interhemispheric surfaces. A semicircle 80 mm in diameter is removed from the iMSP to simulate the SFS. The center of this semicircle is located 20 mm anterior to the center of the basal plane. The posterior ends of SFS and the TI meet at the falcotentorial junction, which is adjacent to the most posterior aspect of the midbrain.
We generate the 3-D FE model of the brain as described above on ANSYS 10.0 ED environment using solid modeling. The X, Y, and Z axes of the Cartesian coordinate system refer to the right-left, anterior-posterior, and cranial-caudal (superior-inferior) anatomical axes, respectively, and the iMSP is represented by the X-Z plane. The completed model of the cerebral hemispheres has only 12 keypoints, 21 lines, 12 areas and 2 volumes. Each cerebral hemisphere has 5 areas in its outer surface, including a convexity surface that is divided into two areas, and a basal surface that is divided into the skull base/tentorium portion and the TI portion, which is subdivided into two areas. The iMSP between the hemispheres is divided into two areas: the falx and the SFS. We have validated this model on CT craniometric data obtained from 50 patients with brain concussion [10].

IMETI 2017
Because the model is only used to simulate deformation of the brain within the cranial cavity caused by intracranial lesions, it consists of only one type of linear elastic material: the brain itself. We omit all extracranial materials, such as the dura, the skull and the scalp, which must be considered in modeling primary brain injuries. Under increased ICP, the cerebrospinal fluid (CSF) spaces such as the ventricles and the cisterns are compressed, so we also neglect all CSF spaces in the model. The specific material properties of the brain are still controversial. From the literature, the Young's modulus (YM) of the brain is taken as 10000 Pa [11] and the Poisson's ratio (PR) is taken as 0.45 [12].
We choose the SOLID187 element, a higher order 3-D, 10-node tetrahedral element, to construct our model. This element has a quadratic displacement behavior and is well suited to modeling irregular meshes [13]. The FE model is meshed automatically at "SmartSize" level 5 using the built-in tool of ANSYS, generating 1391 nodes and 776 elements.

Experimental setup
We perform a single load step, static structural analysis of our brain model simulating a left-sided subdural hematoma (SDH) compressing the cerebral hemispheres. We apply uniform compressive pressure loads on the convexity surfaces of the left cerebral hemisphere, simulating increased ICP caused by the SDH. Since the falx has a much higher YM than the brain [14], we consider it immobile in all degree of freedoms (DOFs). Displacements of the convexity surface of right cerebral hemispheres are also constrained to 0 in all DOFs, simulating the rigid skull covering the brain.
On the basal surface of the brain model, only the TI portion is allowed to move in all directions. The remaining basal surface of the brain lies on the bony skull base anteriorly and the tentorium cerebelli posteriorly. The former structure is fixed to the bone and cannot move. Similar to the falx, the tentorium is also a tough structure and is also considered immobile in our model. As a result, parts of the basal surface of the brain other than the TI are only allowed to move in X-and Y-directions and their displacements are constrained to 0 in Z-direction. The number of surface load conditions is 66 and the number of constraints is 974.
We evaluate the deformation of our model in two sets of experiments. In the first one, we apply different pressure loads ranged from 5 mmHg (667Pa) to 80 mmHg (10664Pa) using material properties derived from the literature. In the second one, we apply a 20-mmHg pressure load, slightly higher than the upper normal limit of ICP [15], to our model having different PR values ranged from 0.25 to 0.49. Different values of YM are also tested. We choose the large displacement static analysis mode to solve our problem. Using the postprocessor of the ANSYS system, the nodal solutions for the displacements are listed and analyzed. Then, the deformed shapes are plotted for evaluation. The displacement over nodes located at the basal circumference (BCC, 40 nodes), the basal midline (BML, 25 nodes), and the horizontal line across the tentorial incisura (HTI, 9 nodes) are recorded. These three groups of nodes are employed as "local approximates" of the maximal displacements of the entire model. Their displacement vector sums as well as the X-, Y-and Z-components are compared to the corresponding global maxima to see if the latter can be estimated.

Brain deformation under different pressure loads using material properties derived from the literature
We apply different pressure ranged from 5 mmHg (667Pa) to 80 mmHg (10664Pa) to the model having a YM of 10000 and a PR of 0.45. Brain deformation after applying a 20-mmHg pressure load is shown in the upper middle part of Fig. 1, using a cutaway view from the anterior-inferior aspect using coronal (Y-Z) cutting planes 20 mm anterior to the center of the basal planes through the center of the subfalcine space. A large proportion of the compressed hemisphere, extending from the outer surface all the way to the SFS, has large displacement, with the subfalcine region where clinicians measure the MLS having the largest.  The magnitude of brain displacement versus pressure load is summarized in Fig. 2 largest total displacements (USUMmax), the vector sums of UXmax, UYmax and UZmax, range from 1.32 mm to 20.05 mm. UXmax along the basal midline of the model (UXmax_BML), which is similar to the MLS, range from 1.19 mm to 16.86 mm. The relationships between these five parameters and the pressure load are almost linear. Nevertheless, UZmax seems to rise slightly slower as the pressure load increases. The two largest parameters, UXmax and USUMmax, increase by 2.5 mm with every 10 mmHg increase in pressure. A 5-mm largest displacement in our brain model (4.97 mm for UXmax and 5.13 mm for USUMmax) occurs with a pressure load of 20 mmHg, which is the treatment threshold recommended in the clinical practice guideline [16]. The difference between UXmax and UXmax_BML ranges from 0.08 mm to 3.00 mm, being less than 1 mm with pressure loads under 40 mmHg.
The accuracy of the local approximates, which are easily visualized on patient images, versus pressure load is summarized in Fig. 3. UYmax along the basal circumference (UYmax_BCC) accurately describes the global UYmax, having ratios larger than 0.99 under all conditions. UZmax along the horizontal line across the tentorial incisura (UZmax_HTI) also accurately describes the global UZmax under all conditions. All UZmax_HTI/UZmax ratios are larger than 0.97. The ratio between MLS (UXmax_BML) and UXmax is 0.93 with a pressure load of 5 mmHg, and then decreases to 0.91 with a pressure load of 40 mmHg. This ratio drops further with even larger pressure loads. The largest X displacement along the BCC (UXmax_BCC) is about 0.89 times UXmax with pressure loads larger than 20 mmHg.

Brain deformation with different material properties
Displacement of our brain model with different PRs after applying a 20-mmHg pressure load is shown in Fig. 1, using cutaway views through the center of the subfalcine space. Values of 0.30, 0.42, 0.45 and 0.49 are used as examples of very low, low, "normal", and high PRs. With a very low PR, the region with the largest displacement in the brain model only occurs at the anterolateral aspect of the BCC. Away from this region, the displacement decreases as the distance increases. With "nearly normal" PRs, a large proportion of the compressed hemisphere has the largest displacement. The region extended from the outer surface to the SFS, showing significant SFH compared to brain shift at the BCC. However, the region with the largest displacement moves abruptly from brain surface to the subfalcine region with PRs larger than 0.434. With a very large PR, the area with the largest displacement is localized inside the TI, showing TTH as the predominant mode of brain shift.
On the other hand, displacement of our brain model is inversely related to the YM given the same PR and the same pressure load. The relative brain deformation is not affected by the YM, so the contour plots with different YMs are not shown.
The magnitude of brain displacement versus different PRs is summarized in Fig. 4. UXmax in our brain model ranges from 14.88 mm with a PR of 0.25 to 3.21 mm with a PR of 0.49. Likewise, UYmax ranges from 8.16 mm to 1.02 mm. On the other hand, UZmax decreases from 7.82 mm to the minimal value of 2.40 mm with a PR of 0.44, then increases gradually to 2.87 mm with larger PRs. USUMmax decreases from 14.91 mm with a PR of 0.49 to 5.64 mm with a PR of 0.43, then slowly decreases to 3.81 mm with a PR of 0.49. The slower decrease in USUMmax is attributed to gradual increase in UZmax. The values and accuracy of the local approximates versus PR is summarized in Fig. 5. In the left part, the absolute values of UXmax, UXmax_BCC and UXmax_BML are plotted against PR. The relationship between UXmax_BCC and PR is linear, while the ratio between UXmax_BML and UXmax_BCC appears to be exponential. UXmax is always close to the larger values among UXmax_BML and UXmax_BCC, with the largest absolute difference of only 0.35 mm. The percentage error is about 5% with larger PRs, when UXmax is approximated by UXmax_BML. In the right part of Fig. 5, the accuracy of the local approximates versus PR is plotted. UYmax_BCC accurately approximate the global UYmax with PRs smaller than 0.48. The ratios

IMETI 2017
between UYmax_BCC and UYmax are 0.95 and 0.67 with PRs of 0.48 and 0.49, respectively. The ratios between UZmax_HTI and UZmax are larger than 0.97 with PRs larger than 0.434, revealing the central position of the UZmax nodes. UZmax_HTI becomes a poor approximate when the node having UZmax is at the surface of the compressed cerebral hemisphere. Likewise, the MLS (UXmax_BML) can approximate UXmax better with larger PRs, when the node having UXmax is near the TI. On the other hand, UXmax_BCC can approximate the UXmax with PRs smaller than 0.434. The ratio between UXmax_BCC and UXmax decreased to 0.96 with a PR of 0.44, then decreased further and is surpassed by UXmax_BML/UXmax with larger PRs.

Discussion
We have proposed a simplified FE model that can be applied to simulating brain shift caused by an SDH. Our model is constructed with a rather small number of nodes and elements, but is able to simulate brain displacement and herniation in real patients. The pressure loads we tested covered almost all ICP conditions encountered clinically and the resulting patterns of brain shift appear similar to those appeared on images in patients with acute or chronic SDH, as shown in Fig. 1.

Meaning of each parameter
In the first part of our experiment, our model with "standard" biomechanical parameters is subject to different pressure loads. Pressure gradient across different compartments of the intracranial space does exist under certain conditions such as posterior fossa tumor [17], but the ICP in patients with intracranial hematoma usually equalize within hours. Clinically, the pressure is almost always considered uniform within the dura, and it is acceptable to place the monitor anywhere. To date, there is no specific recommendation regarding the best location of the ICP monitor in patients with severe traumatic brain injury [18]. We also follow this approach by assuming that there is only one equalized pressure value within the intracranial space.
In the second part of our experiment, models with different PRs are subject to a pressure load slightly higher than normal ICP. We have tested these models with different pressure loads, and the relationship between the maximal displacements and the pressure appear to be linear in the range between 5 and 80 mmHg, obviating the need for testing with multiple load conditions. With all PRs, the UXmax node at the outer surface of the compressed cerebral hemisphere always locates at the BCC. Therefore, UXmax_BCC can be regarded as the point with maximal SDH thickness. Brain models with low PR have MLS is much smaller than the UXmax_BCC. They behave like atrophic brains compressed by chronic SDH, a condition commonly seen in aged patients. On the other hand, models with very high PR resembles the behavior of tight or "angry" brains, where small volumes of SDH, represented by small UXmax_BCCs, causes larger brain deformation at its center and squeezes it through the TI, as noted in the autopsy and image findings with TTH [15]. From the regression formula derived from Fig. 5, the thickness of the SDH ranges between 0.15 mm with a PR of 0.50 and 29 mm with a PR of 0, under an ICP of 20 mmHg, concordant with clinical observations. Although many complex constitutive models had been proposed previously, Wittek et al. had pointed out that the choice of the brain tissue constitutive model, when used with an appropriate finite deformation solution, does not affect the accuracy of computed displacements, and therefore a simple linear elastic model for the brain tissue is sufficient [19]. Our results support their opinion. In our experiments, we are able to use a single global variable, the PR of the simple linear elastic model, to model the "tightness" of the whole supratentorial brain to produce realistic deformations. There is no other biomechanical model in the literature that can simulate both acute and chronic SDH without any change in its geometry to accommodate the loss of brain volume.
The "slower" UZmax change to pressure load is related to the relatively coarse meshing of our model. In this paper, the number of nodes/elements is 1391/776. Using a finer mesh consisted of 12450 modes and 8263 elements generated with a full version of ANSYS, UZmax changes almost linearly with pressure load up to 70 mmHg.

Clinical implications: subfalcine and transtentorial herniations
We have carried out several studies on automatic recognition of the MLS, the evidence of SFH on brain CT images [4][5][6]. Using visual inspection, we use a quadratic curve to model the deformed midline at the foramen of Monro, near the base of the cerebrum. In the simple FE model proposed in this paper, the shapes of all deformed basal midlines under different test conditions assume the shape of quadratic curves as we observed on patient images. This model may provide a biomechanical basis for modeling the deformed midline and quantifying the MLS.
We set the value of PR as 0.45 for the simple linear elastic model in this study [12]. To verify this value, we measured the thickness of the SDH and the corresponding MLS on initial CT images in 12 surgically-treated acute SDH patients younger than 50 years of age, collected from our previous studies [5,20]. All patients showed signs of increased ICP. Regression analysis showed a linear relationship with a Pearson's correlation coefficient (R2) of 0.956. The MLS is about 0.96 times the thickness of SDH (95% confidence interval: 0.83 to 1.10). This ratio between the MLS and SDH thickness is very close to our FE model with a YM of 10000 and a PR of 0.45, which have MLS/SDH ratio between 0.95 and 1.03 with pressure loads larger than 30 mmHg, when surgical decompression is necessary. To fully validate our model, large-scale studies incorporating CT images and corresponding initial ICP data measured intra-operatively are required.
In a model with PR close to 0.45, the region with the largest brain shift extends from the anterolateral region of the affected cerebral hemisphere toward the subfalcine region. This region is decompressed by current practice of unilateral front-parieto-temporal craniectomy for acute SDH with or without brain contusion [16,21]. From our FE experiments, we also show that MLS is a good estimator of USUMmax in models with medium to large PRs (Fig. 5, right). In addition, the MLS with a 20 mmHg pressure load never exceeds 7 mm even with very low PRs (Fig. 5, left), showing the robustness of MLS as a sign of increased ICP, especially compared to the shift at the brain surface, or the size of the mass.
In addition to SFH, our model also demonstrates TTH, shift of the brain region adjacent to the TI. Compared to the SFS, the TI is relatively small in area. As a result, qualitative and quantitative description of TTH is much more difficult than describing the MLS [2]. Change in brain stem width had been proposed as a possible sign of compression [22]. In our experiments, the UXmax_HTI, behaves similar to the MLS, i.e. UXmax_BML with different ICPs and PRs.
On the other hand, UZmax does increase with PRs larger than 0.44 (Fig. 4). Our model is more consistent with the classic observation, which describes TTH as downward movement of the brain stem. However, the magnitude of UZmax_HTI is rather small. Even with very large ICP and very large PR, the maximal UZmax_HTI never exceeds 7 mm, which is far smaller than the largest possible MLS. The result is concordant to the case report by Ropper in 1993 [23]. Using MRI, the author is able to demonstrate the 2-mm downward displacement and the 13-mm horizontal displacement of the upper brain stem.

Limitations of this study
Our study has several limitations. The simple FE model is an "average" across subjects and cannot be used to provide the detailed deformation of each part of the brain on specific patient without modification. To simplify our model, the internal structures of the brain are omitted. As a result, detailed point-to-point correspondence, even to the nonlinearly transformed anatomical atlas, is not possible. We assume the PR of the model to be constant, while the brain may change its PR over time, either locally or globally, over time when it is compressed by the hematoma.