Failure Criteria for Prediction of Clinically Relevant Damage of Annulus Fibrosus

Internal disc damage and annular tears are some of the causes of back pain. It is difficult to observe location of internal damage and its distribution and propagation in cadaveric experiments or via imaging in clinical scenarios. Finite element analysis is useful for understanding the effects of different loading conditions on the location of initiation and propagation of tears in the annulus. In this study, a hyperelastic anisotropic material model in conjunction with biaxial properties was used for modelling annulus fibrosus (AF) material under a variety of loading conditions. The loading conditions considered included compression, flexion, extension, and their combinations. Maximum normal stress (MNS) and Tsai-Wu criteria were used to evaluate the damage. Strength anisotropy and asymmetry were considered in both criteria. The damage predictions based on Tsai-Wu criterion were in better agreement with available experimental studies and clinical observations.


Introduction
Internal disc damage (IDD) is a cause of chronic pain that 80 percent of people experience sometime during their life. It decreases quality of life and limits physical activities [1]. Slightly degenerated discs in the lower lumbar intervertebral disc (IVD) are more susceptible to disc prolapse, associated with back pain [2]. Formation of tears in the annulus fibrosus (AF) is recognized as the first sign of the disc degeneration process [3,4]. The tendency of the disc to herniate posteriorly and posterolaterally depends on many factors, including the loading conditions. Theoretical and experimental attempts have been made to address damage formation in the annulus. However, the relationship between initiation location and propagation of tears and loading conditions is not entirely explained.
With respect to the clinical studies and MRI images, there are three categories of annular tears in AF [5]: radial, circumferential (concentric), and rim (peripheral) lesions (Figure 1). Radial tears have been observed in the posterior and postero-lateral regions of the lower lumbar spine [6]. Concentric tears were distributed equally in the anterior and posterior regions [3]. Discrete peripheral tears or rim lesions [7] were seen more frequently in the anterior region at the junction of the endplate and disc [3,4]. The combination of these tears might also be seen in a disc, depending on the age, genetic background, and load history.
Physiological loading produces gradual prolapse. It starts with distortion in the lamellae of AF and forms fissures which cause protrusion of the nucleus polposus (NP) [8]. Experimental attempts involving human cadavers and animal models were made to investigate the relationship between annular tears and loading conditions [8]. Specimens were subjected to cyclic axial compressive loads ranging from 37 to 80 percent of their failure This article is availabel in: http://spine.imedpub.com/ Spine Research ISSN 2471-8173 load limits and discs damage was reported in 2 of 11 specimens tested [9]. Vertebral or endplates damage were reported in the remaining 9 of 11 specimens [9]. In another study, cyclic compressive loads were applied to a motion segment with a predefined flexed position in which distortion of the posterior lateral area of the disc was observed [10]. Flexion-extension in addition to compression at three load levels of 260 N, 867 N, and 1472 N produced posterior and postero-lateral disc herniation [11]. Combination of flexion, axial rotation and compression produced annular tears and disc prolapse, which were peripheral at origin [12].
Conducting experiments to study damage initiation and progression are valuable to understand the injury mechanism of AF. However, it is difficult to identify the damage and its location and propagation under different and combined or multi-axial loading conditions. Moreover, these experiments may have error caused by uncontrolled variations of experimental conditions. Finite element (FE) approach is an effective tool which has been used for understanding spine kinematics in response different loading conditions [13][14][15][16][17]. Damage of the disc and stress distribution on the annulus fibrosus due to the different loads can also be studied using FE technique.
Damage and failure criteria are used in conjunction with FE analysis by several investigators. These criteria are based on strain, stress, or energy state within the AF. The maximum tensile strain criterion of the collagen fibers was used for failure analyses in AF under single and combined loads [18,19]. In another study intradiscal pressure in the nucleus, shear strain between annulus and adjacent endplates, and fiber strain in the annulus under combined loading scenarios were considered for predicting damage location [20]. The failure stress and strain data for soft tissues are quite inconsistent in the literature, depending on the strain rate used in performing the experiment. Three envelopes of yield and failure for collagen fibers were produced based on strain, stress, and strain energy density (SED) criteria using a nonlinear viscoelastic FE model of fibers [21]. The authors concluded that the critical value of stress did not change for a wide range of strain rates when SED criterion was used. However, ultimate stress and strain were different for a wide range of strain rates in various regions of the functional spinal unit. It was concluded that SED is a more practical failure criterion for AF region [21].
A disc model with stage (I) degeneration was used to explain delamination in the AF lamellae based on interlaminar shear stresses [22]. Maximum principal strain distribution in AF region was used to evaluate damage in a model based on hyperelastic anisotropic constitutive behaviour of AF [23]. Another study developed models of a healthy and a degenerated disc, to predict delamination using principal strain distribution [24]. In a recent study, damage initiation and progression in AF region under various cyclic loading conditions was studied using a 3D poroelastic FE model. The damage parameter used was the maximum principal tensile stress criterion in the annulus ground substance, along with linear damage accumulation for failure propagation [25]. Poro-elastic FE was developed to predict failure progression in a lumbar disc under heavy lifting cyclic loads [25]. The effects of fluid flow and strain dependent permeability were included in the analyses. Failure was defined when von Mises stress exceeded 2 MPa in the annulus matrix, 20 MPa in the nucleus and 4 MPa in the endplate [26].
Previously discussed studies either used a certain component of stress/strain such as fiber strain to quantify damage in the annulus under multi-axial stress/strain states, or combined different components of stress/strain in the form of von Mises distortion energy or principal stress/strain. However, these criteria have not been modified to account for the anisotropy observed in the AF material properties, especially its strength. Moreover, there is a limited knowledge about propagation of damage under combined loading conditions out of these studies. The Tsai-Wu criterion, a generalized form of von Mises distortion energy, has been used for anisotropic biomaterials such as tooth tissue [27] and vertebral trabecular bone [28]. Therefore, it is the goal of this work to apply and extend the use of this criterion to anisotropic soft tissues, such as AF.
Biaxial material properties of AF utilized in a hyperelastic anisotropic material model were shown to provide accurate results in terms of the biomechanical behaviour and stress distribution in the first part of this study [29]. Therefore, biaxial material properties of AF were used in present study. The aim of was to use two appropriate failure criteria for anisotropic materials, Tsai-Wu and maximum normal stress, to evaluate damage initiation and propagation under different loading conditions. Tsai-Wu failure criterion considers stresses and strengths of the material in different directions and results in a scalar value. Results from the Tsai-Wu criterion are then compared to maximum normal stress criterion to evaluate damage location and distribution in AF and its clinical relevance in predicting annular cracks.

Finite element model
A previously validated 3D FE model of ligamentous functional spinal unit (L4-L5) [30] was modified for the analyses. The inferior surface of the L5 vertebra was rigidly fixed [29]. A hyperelastic anisotropic material model available in ABAQUS/Standard TM version 6.11 (Simulia, Inc., Providence, RI, USA) referred to as Holzapfel-Gasser-Ogden (HGO) [31,32] model, was used as the material model for AF. Biaxial material properties were used as material constants of AF for damage analyses. Properties of other elements of FE model such as bony elements, ligaments, and NP were kept the same as in the previous model [29]. Equation (1) shows the governing strain energy relation in the AF region: where C 10 , K 1 and K 2 are material parameters and the matrix is considered as incompressible material. Detailed description of this material model and values of constants for both uniaxial and biaxial conditions are reported in [29].

Strength properties in different directions
Strength varies as a function of direction in AF. Therefore, for conducting damage studies, stresses are compared to their respective strength limit of the annulus in the corresponding direction. Tensile experiments in the uniaxial state were performed in axial [33], circumferential [34], and radial [35] directions on multiple layers of AF lamellae and strengths of this tissue in these orientations were reported. Different anatomical regions of AF (posterior inner/outer and anterior inner/outer) represent different strengths and fracture strains in each direction as listed in Table 1. Mean experimental values of strengths in the axial and radial directions and mid-range value of circumferential direction were used in this study. Posterior and anterior regions were separated by a plane connecting the left and the right lateral regions. Three inner layers of AF close to the NP were considered as the inner region, and the remaining layers considered as the outer region, (Figure 2).
Literature on the compressive strength data of AF is very sparse. An experimental procedure was proposed [36] to study the mechanical response of this tissue to both initial and repeated compressive axial loads. Samples were prepared such that no continuous collagen fibers were coupled to the two endplates of the specimens, therefore, collagen fibers were not actively bearing load [36]. In the HGO material model used in this study, fibers cannot sustain compressive stress. Therefore, the compressive strength of AF was assumed as compressive strength of the ground substance from the experiments [36]. Annulus ground substance which is considered as an isotropic material [15][16][17][18][19][20][21][22], has the same strength in all directions. The compressive stressstrain curves up to failure of the annulus for different regions were obtained under monotonic and cyclic loading and similar failure compressive strengths for different regions of about 2 MPa were obtained under both types of loadings [36].
There are no published experimental results regarding the shear strength of AF. An analytical method was used to calculate shear strength from uniaxial tensile experiments on a single layer of AF [37]. Experiments were performed in longitudinal, transverse, radial, and circumferential directions. Based on Tsai-Hill failure theory for a symmetric composite, the shear strength was calculated to be 0.24 MPa. This value was used in the present study.

Failure criteria
Two commonly used criteria to evaluate damage under multiaxial stress states are maximum normal stress and distortion energy criteria. However, more general forms of these criteria are needed for anisotropic materials such as AF. These are described in this section.

Maximum normal stress criterion
There are three planes of principal stresses on which there is no shear stress. To find tensile strength of an anisotropic material in a particular direction (n), the following stress transformation is used: (2) where, s are the material tensile strengths in three orthogonal directions. Axial, circumferential, and radial strengths of AF were considered as characteristic orientations (Figure 3a). Regional variation of the tensile strength in the anterior versus posterior and inner versus outer were considered based on the values listed in Table 1. Angles θ and ϕ specify the principal directions for each element, as illustrated in Figure 3b. Strength asymmetry of the material was considered and compressive strength in all directions was considered to be 2 MPa, as explained earlier.
To use maximum normal stress criterion the following steps were performed: a. Stress components from FE simulations for each element of AF were used in a Matlab code to find principal stresses and directions for each element.  --*Mid-range values were calculated in this study based on the range provided in the reference. Table 1 Tensile strength and fracture strain in axial, circumferential and radial directions of the different regions of AF.

c. Absolute values of maximum and minimum principal stresses
were divided by the corresponding tensile or compressive strength of the material in the same principal direction. The maximum value of the two ratios for each element was considered as damage ratio. Damage distribution in AF for each specific loading condition was then obtained.

Tsai-wu criterion
The Tsai-Wu quadratic criterion can account for strength asymmetry and anisotropy, as well as interaction between strengths under multi-axial loading conditions. This failure criterion is a generalized form of von Mises criterion. The basic assumption of Tsai-Wu criterion is that there is a failure surface in the scalar form which can be found by the following equation [39]: where F i and F ij are the second and fourth rank strength tensors, respectively. The Tsai-Wu equation may be expanded and written in term of damage (d) for no failure, in the following form: (4) where s are normal stresses in three perpendicular directions and s are shear stresses. Stresses with respect to axial, circumferential, and radial directions for each element were used in this equation. The strength dependent F values in the equation were computed based on the assumed strengths accounting for the difference in tensile and compressive properties of different directions. Tensile properties for various anatomical regions of AF were considered, as listed in Table 1 and compressive strength was considered to be 2 MPa, as described earlier. Values of coefficients F's and their corresponding equations are listed in Table 2. In this table, X, Y, and Z are tensile strengths in axial, circumferential, and radial directions, respectively, and X', Y', and Z' are the corresponding compressive strengths. Since compressive strength was assumed to be 2 MPa in all directions as previously discussed, X' = Y' = Z'. Also, Q, R, and S are shear strengths on circumferential-radial, axial-radial, and axial-circumferential planes, respectively, and Q', R', and S' are the corresponding shear strengths in the opposite direction. All of the shear strength values were considered to be 0.24 MPa, as discussed earlier. The contribution of shear stresses on axial-radial and circumferential-radial planes in the Tsai-Wu damage criterion is small, in comparison to axial-circumferential shear component. This is because of the small ratio of these stresses when divided by the corresponding shear strength of the material.
To use the Tsai-Wu failure criterion for AF, the following steps were executed: a. Stress components in axial, radial, and circumferential directions for each element were obtained for each load condition. Analyses were done for different transverse surfaces in the height direction (inferior, superior, and mid-height). These surfaces were chosen based on the available planes with experimental results or the plane where maximum damage occurred.

Loading Conditions
It has been shown that follower loads can model compression due to muscle forces in the human spine, similar to the physiological conditions [40]. Therefore, all the compressive forces were simulated using follower loads concept in this study. Bending moments were applied to a reference point which was connected to all of the nodes of the superior surface of L4 by kinematic coupling. To assess damage location, various simple and combined loading conditions were applied to the functional spinal unit FE model (FSU), as listed in Table 3. Simple loads included compression, extension, or flexion. Combined loads included combination of compression with flexion or extension.
A compressive load of 500 N was considered for the analysis, as suggested for modelling physiological loading condition of the L3-L4 level segment during standing posture [41] (load case 1, ( Table  3). A compressive load of 4000 N has produced complete disc failure during cyclic loading in experiments [10]. Therefore, this load level was also considered to study propagation of damage in AF (load case 2, Table 3). Pure flexion and pure extension within the physiological range of intact motion segment, corresponding to 7.5 N.m moment were also applied to the model (load cases 3 and 4, Table 3).
For the combined loading scenarios, a 500 N compressive load and 7.5 N.m bending moments in flexion and also in the extension were applied to simulate loading close to the normal physiological condition (load cases 5 and 6, Table 3). It has been noted that during lifting activities compressive load that is applied to the lumbar spine by back muscles can reach six to eight times the weight which is lifted [10]. Therefore, a 4000 N compression load was applied to a hyper-flexed motion (8° flexed) segment to model heavy lifting similar to the available experiment studies [2][3][4][5][6][7][8][9][10] (load case 7, Table 3). It has also been shown that human motion segments during actual flexion is oriented at the x-axis minus 30° (based on the terminology of Panjabi et al. [42]) [43]. This loading condition with a combination of 1400 N compression Figure 3 (a) Arrangement of fibers in the annulus fibrosus with respect to the axial, radial and circumferential coordinate [38], and (b) difinition of angles in Equation (2) for the strength transformation in principal directions.

Spine Research ISSN 2471-8173
(moderate lifting) and 7° flexion in this orientation on the L4-L5 disc was also considered (load case 8, Table 3). Except for results of natural flexion loading condition which was not symmetric, other loading conditions were symmetric with respect to the mid sagittal plane and, therefore, half of the disc is shown in the figures.
In another scenario, combination of 500 N compression load and 7.5 N.m flexion moment (5° rotation) in one model and 500 N compression load and 7.5 N.m extension moment (2° rotation) in another model were considered as initial conditions. Then in both models first the magnitude of compression force was increased to 2000 N (this is considered to be a high level of compression force) and the increase in damage values were calculated and compared to the initial model. Second, the compression load was kept at 500 N and the magnitude of flexion or extension in initial models was increased to the level that produced 8° or 4° rotation (considered to be a high level of rotational angle [40]), respectively, and the increase in damage value were compared with the initial models.

Damage distribution and initiation location
Based on the maximum normal stress criterion (MNS), under a compressive force of 500 N, damage initiates at the anterior outer region of the inferior surface (Figure 4b). Inferior surface was chosen to evaluate the results because available experiments for comparison are visible for this surface. Under a higher level of compressive force at 4000 N, damage propagates peripherally along the lateral region and extends radially to the next lamellae in the lateral site (Figure 4d)   -Maximum shear strain was at Posterior inner [20] -Maximum Tensile strain was at Postero-lateral outer [20] 4 Pure extension   When the motion segment is subjected to 1400 N compression and 7° actual flexion angle, based on Tsai-Wu criterion the maximum damage is at the postero-lateral inner part of the annulus and close to the inferior surface, (Figure 9a). According to MNS criterion, however, the highest damage is at the anterior region of this surface, which propagates more peripherally and tends to form a rim tear, (Figure 9b). Damage distributions at the inferior surface of the annulus are compared with the available experimental result for a cut through specimen in Figure 9c. Since the damage distribution was similar in the height direction, inferior surfaces which have The highest damage value are shown. As can be seen from this figure, Tsai-Wu prediction better matches with the experimental result. Therefore, the two criteria predict different damage initiation and propagation and MNS criterion predicts a higher degree of damage. Increasing of damage due to increasing magnitude of pure compression load, pure flexion and extension moments on FSU based on the MNS and Tsai-Wu criteria is shown in Figure 11. This figure indicates that according to both criteria flexion produces a higher degree of damage than extension. However, this difference is larger for Tsai-Wu criterion. It can also be seen that while compression load damage increases nonlinearly with increasing load, for flexion and extension the trends are almost linear. Under pure compressive load, Tsai-Wu criterion predicts higher damage value than MNS criterion. However, with increasing the load level, damage values of the two criteria become closer.

Comparison of loading cases and effect of increasing load
Effects of compression and bending moment (flexion or extension) in combined loading conditions (load cases 5 and 6) are shown in  locations, respectively ( Figure 12). However, by increasing flexion and extension moments to a rotation of 8° and 4° in loading cases 5 and 6, respectively, increase in damage values are similar at 25% (Figure 12). Based on either criterion, it can be concluded that higher compressive forces produce more damage than higher bending moments. Also, increasing compression in extended FSU (load case 6) produces higher damage value than flexed FSU (load case 5) using either criteria.

Discussion
During different complex loading scenarios of spine various tissues can be at risk of injury. Damage in the disc can occur due to high level of loading or because of relatively low loads that are repeated many times. Herein, we investigated the risk of failure for different regions of annulus fibrosus under diverse static loading conditions to identify the location of damage initiation and its distribution by applying two damage criteria. Findings of this study correspond to and are validated by the damage location observed in experimental and clinical studies available in the literature.
During daily activities compressive forces, including body weight and muscle forces, are applied to the spine. These compressive loads increase to higher levels during lifting exertions. The Experimental results indicate that under high magnitude of compression endplates damage and rim annular tears are more likely to occur in the intervertebral disc [44]. Degeneration and damage of disc due to low level of compression produces posterior and postero-lateral radial damage [7] as predicted by Tsai-Wu criterion.
In-vitro experiments show that when the spine is in a flexed position the nucleus moves posteriorly and when the spine is in the extended position it moves in the opposite direction [48], which can produce higher stresses in those regions. Based on both damage criteria, under pure flexion posterior region of the disc and under pure extension anterior regions were critical locations. However, disc failure is unlikely to occur by bending moment in the absence of compressive loads and under this loading condition mostly ligaments are at risk of failure [25] which was not the focus of this study. In addition, body weight is always applied to the spine motion segments in-vivo conditions. Under combined bending moments of 7.5 N.m and 500 N compression, critical location moves to the posterior and postero-lateral inner region of the annulus according to Tsai-Wu criterion, which also matches with other relevant FE results [25]. However, anterior outer and posterior outer damage locations under combined extension-compression and flexion-compression can be explained by clinically observed rim tears [7]. Based on this analysis, it can be concluded that during lifting objects it is better to keep a slightly flexed position. This can be explained by a lower intradiscal pressure during this posture than in upright position.
High compressive load such as 4000 N can be applied to a fully flexed motion segment, having 8° forward rotation while landing [49]. Damage analysis of this study by using both criteria predicted damage at the postero-lateral region of the inferior surface of AF in the form of radial tear, supported by in-vitro experiments [3] and clinical observations [50]. Based on the analyses, compressive force has the largest effect on the developing damage in AF. This situation happens mostly while lifting heavy items, which is one of the most common loading conditions leading to 50 to 75 percent of the low back pain [51]. Weight of an object is the most important factor for controlling damage during lifting.
Combination of flexion and compression is a common movement in real working life as documented by Fathallah et al. [52]. During these tasks, there is a coupled movement in the lateral direction [43], which was considered in this study by 30° orientation in the flexion direction with respect to the anterior-posterior axis. The initiation and propagation of herniation were studied using nucleus tracking in-vitro experiment [43], as shown in Figure 9c. These results indicate that herniation mostly initiated and propagated posterolaterally [2]. This is predicted by Tsai-Wu criterion for natural flexion combined with a moderate compressive load. Although MNS criterion for combined forward flexion and compression predicted clinically relevant results, its prediction was different from available experimental result [43] for natural flexion. Based on the results of this study and according to Tsai-Wu damage criterion, except for extension, accumulation of damage is at the posterior region of the disc and close to the inferior endplate. High intensity zones (HIZ) in MRI images that are considered as a sign for diagnosing internal disc disruption [53], were observed posteriorly and posterolaterally (Figure 13) in people with low back pain [54,55]. MNS criterion does not appear to accurately predict failure due to the accumulation of damage under daily physiological loading conditions. In term of propagation of annular tears also Tsai-Wu seems to present reasonable results. Radial tears were observed in compression, combination of compression and flexion, and combination compression and flexion in 30° orientation (combined flexion and lateral bending). Concentric tears or delamination of the annular was not observed in the MRI images [55], whereas these are known as part of the damage during degeneration process. In all of the loading conditions, damage according to both criteria propagated in the circumferential direction due to the high shear stress [23] in the middle layers.
Several assumptions and limitations of the analysis in this study should be noted. Water content and porosity of the disc were not included in the analysis, although they can have an impact on the damage evaluation of soft tissues. However, for the duration of static loads considered in this study these effects are thought to be small. Further analysis needs to be conducted to evaluate these effects. The healing ability of in-vivo organs is also difficult to understand with current knowledge and it was not incorporated in this study. The annulus region was considered as a homogeneous material in term of fiber orientation and density. It has been shown that the orientation of fibers changes from the inner to the outer part of AF [56]. However, this assumption should not change the results obtained significantly, since different strength for different anatomical regions were considered. For more accurate analyses, it is necessary to consider this heterogeneous property of annulus. Viscoelasticity is another characteristic of soft tissues [57,58] which was not considered since the focus of this study was on damage initiation location and propagation. Finally, in this study strength of the material was considered under static loading. However, for calculating life of the tissue, although challenging, extension of the results to cyclic loading conditions is necessary.

Summary
Although a number of experimental and numerical studies over the last two decades have significantly improved our knowledge of damage in AF, consensus on a proper mechanical parameter to quantify the damage does not yet exist. In this study, Tsai-Wu and maximum normal stress (MNS) criteria were used to evaluate and quantify damage under a variety of loading conditions including compression, flexion, extension, and their combinations. The FEA results obtained in this study were compared with available experimental and clinical observations from the literature. The two damage criteria generally predicted different damage initiation locations, propagation, and distribution under the different loading conditions. MNS criterion results in a higher damage value for the great majority of the loading conditions considered. Tsai-Wu criterion predictions were in general in better agreement with the available experimental results and clinical observations. Based on both damage criteria, for the same moment magnitude a higher damage is caused in pure flexion than in pure extension, and combination of flexion and low level of compression decreases the maximum damage value compared to pure flexion. Also, higher compressive forces produce more damage than higher bending moments and increasing compression in extended FSU produces a higher damage than a flexed FSU. Based on the analysis presented, accumulation of damage under daily physiological loading conditions is at the posterior region of the disc, in accordance with the Tsai-Wu damage criterion and in agreement with MRI images in people with lower back pain. It also seems that compression plays an important role in appearance of radial tears in combined loading conditions and by applying flexion, damage is more concentrated in the posterolateral region of the AF. In all the loading conditions concentric tears are observed and explained by delaminated AF in almost all of the degenerated discs.