Predictive Factors of Kinematics in Traumatic Brain Injury from Head Impacts Based on Statistical Interpretation
Xianghao Zhan, Yiheng Li, Yuzhe Liu, August G. Domel, Hossein Vahid Alizadeh, Zhou Zhou, Nicholas J. Cecchi, Samuel J. Raymond, Stephen Tiernan, Jesse Ruan, Saeed Barbat, Olivier Gevaert, Michael M. Zeineh, Gerald A. Grant, David B. Camarillo
NNoname manuscript No. (will be inserted by the editor)
Predictive Factors of Kinematics in Traumatic BrainInjury from Head Impacts Based on StatisticalInterpretation
Xianghao Zhan · Yiheng Li · YuzheLiu · August G. Domel · Hossein VahidAlizadeh · Zhou Zhou · Nicholas J.Cecchi · Stephen Tiernan · JesseRuan · Saeed Barbat · Olivier Gevaert · Michael Zeineh · Gerald Grant · DavidB. Camarillo
Received: date / Accepted: date
Abstract
Brain tissue deformation resulting from head impacts is primarilycaused by rotation and can lead to traumatic brain injury. To quantify braininjury risk based on measurements of accelerational forces to the head, variousbrain injury criteria based on different factors of these kinematics have beendeveloped. To better design brain injury criteria, the predictive power of rota-tional kinematics factors, which are different in 1) the derivative order, 2) thedirection and 3) the power of the angular velocity, were analyzed based on dif-ferent datasets including laboratory impacts, American football, mixed martial
Xianghao Zhan, Yuzhe Liu, August G. Domel, Hossein Vahid Alizadeh, Zhou Zhou,Nicholas J. Cecchi, David CamarilloDepartment of Bioengineering, Stanford University, 94305, CA, USA.E-mail: [email protected], [email protected], [email protected],[email protected], [email protected], [email protected], [email protected] Li, Olivier GevaertDepartment of Biomedical Informatics, Stanford University, 94305, CA, USA.E-mail: [email protected], [email protected] TiernanTechnological University Dublin, Dublin, Ireland.E-mail: [email protected] Ruan, Saeed BarbatFord Motor Company, 3001 Miller Rd, Dearborn, MI 48120, USA. E-mail: [email protected],[email protected] ZeinehDepartment of Radiology, Stanford University, 94305, CA, USA. E-mail:[email protected] GrantDepartment of Neurosurgery, Stanford University, 94305, CA, USA. E-mail:[email protected] Zhan, Yiheng Li and Yuzhe Liu contributed equally to this work. a r X i v : . [ phy s i c s . b i o - ph ] F e b Xianghao Zhan et al. arts (MMA), NHTSA automobile crashworthiness tests and NASCAR crashevents. Ordinary least squares regressions were built from kinematics factorsto the 95% maximum principal strain (MPS95), and we compared zero-ordercorrelation coefficients, structure coefficients, commonality analysis, and dom-inance analysis. The angular acceleration, the magnitude and the first powerfactors showed the highest predictive power for the laboratory impacts, Amer-ican football impacts, with few exceptions (angular velocity for MMA andNASCAR impacts). The predictive power of kinematics in three directions(x: posterior-to-anterior, y: left-to-right, z: superior-to-inferior) of kinematicsvaried with different sports and types of head impacts.
Keywords traumatic brain injury · head impact · regression interpretation · commonality analysis · dominance analysis Traumatic brain injury (TBI), a primary cause of death and disability, posesa global health threat and affects over 1.7 million children and adults in theUnited States and over 55 million people worldwide [1]. TBI may arise fromfalls, vehicle accidents, and popular contact sports such as American footballand mixed martial arts (MMA) [2,3]. The prevalence of TBI and the severeconsequences suffered by TBI patients calls for better monitoring of braininjury risks because early detection and intervention are essential for recovery[4]. To evaluate brain injury risk, several brain injury criteria (BIC) have beenproposed based on reduced-order mechanics models [5,6,7,8] and deep learningmodels [9,10]. However, great disparities are noted among the mathematicalforms of the existing kinematic-based BIC in terms of the factors (combi-nations of kinematic features) [11,12,13], such as: 1) the derivative order ofangular velocity, e.g., the Brain Injury Criterion (BrIC) includes the peak val-ues of angular velocity [12] and the Rotational Injury Criterion (RIC) includesthe integral value of the magnitude of angular acceleration [13]; 2) the differentcomponents of kinematics, e.g., BrIC includes the kinematics in three spatialdirections [12] while RIC includes only the magnitude [13]; 3) power of kine-matics, e.g. RIC includes the integral values with a power of 2.5 [13]. Althoughvarious factors have been included in different model designs, the contributionof each of the factors is not yet clearly understood in this field.Brain injury occurs when the inertial force during head movement causeslarge deformation in the brain tissue. Therefore, brain strain, particularly max-imum principal strain, is regarded as a key parameter indicating brain injuryrisk for TBI research [5,6,7,14]. Previous studies have evaluated the BIC bycomparing it with the 95th percentile of maximum principal strain (MPS95)[5,6,7,9,15]. However, the predictive power of the different kinematic factorsin the regression of the MPS95 requires further investigation. Therefore, thepresent study aims to statistically analyze the predictive power of the differentkinematic factors in the regression of brain strain. itle Suppressed Due to Excessive Length 3
In this study, the kinematics from simulated head model impacts, Americanfootball impacts, MMA impacts, car crash impacts, and racing car impactswere used. On each dataset, four statistical interpretation methods were used:zero-order correlation coefficients, structure coefficients, commonality analysis,dominance analysis. The first two methods analyze the predictive power ofindividual kinematic factors and the last two methods analyze the predictivepower of combinations of kinematic factors. The features used in the linearregression were derived from angular velocity, because linear acceleration hasbeen proved to contribute less to brain strain [16], and the ground-truth MPSwas calculated by finite element analysis (FEA).We analyzed the contributions of the kinematic factors based on: 1) thederivative order (angular velocity, angular acceleration, angular jerk); 2) thecomponents in three directions in the anatomical reference frame and themagnitude; and 3) the power (square root, quadratic, cubic, etc.).
Xianghao Zhan et al. ω (0) ( t ), ω (1) ( t ), ω (2) ( t ). The numerical difference was calculated by a5-point stencil derivative equation. This analysis was limited to these threederivative orders because they have definite physical meanings, i.e., the angu-lar velocity, the angular acceleration, and the angular jerk. Furthermore, thenumerical differentiation generally introduces noise, which makes the high-order derivative of the angular velocity inaccurate.2) The components of rotational parameters were derived from the 3Dkinematics in the anatomical reference frame (x: posterior-to-anterior, y: left-to-right, z: superior-to-inferior). The magnitude of the kinematics was alsocalculated and used as the fourth component. It has been shown that impactsfrom the side and rear directions might be more likely to lead to obviousperformance decrement (OPD) in football players [26], which suggests thatthe specific directions of 3D kinematics should be included in the regressionof brain strain.3) To capture the non-linear relationship between the features and MPS95,powers of 1 , . , , , , , itle Suppressed Due to Excessive Length 5 tures and the corresponding MPS95 of each impact were visualized in the ex-ample heatmap (Fig. 1) that described dataset NHTSA. (Heatmaps of otherdatasets were showed in Supplementary Figs. 1-6.)2.3 Regression model of MPS95 and model interpretationOrdinary Least Squares (OLS) was used to build linear regression models toquantify the contributions of kinematic factors (predictors) to the predictionof MPS95 (outcome). The OLS model captures the most direct and explicitpredictive power of a factor in terms of explained variance, and therefore itwas used to investigate the predictive power of factors for MPS95.To analyze the relative predictive power of kinematic factors (differentderivative order/component and magnitude/power) in linear regression mod-els, four statistical interpretation methods were used: zero-order correlationcoefficients [28], structure coefficients [29], commonality analysis [30] and dom-inance analysis [31]. Considering the possible collinearity among features se-lected, the methods were chosen because they are interpretation methods with-out the requirement of independence of the predictors. Besides the predictivepower analysis of the direct effect of each factor (zero-order correlation coef-ficients and structure coefficients), the common effects of factors and domi-nance relationships of factors were also studied by commonality analysis [30]and dominance analysis [31], respectively.Zero-order correlation coefficients measure the direct effect by the Pearsoncorrelation between each predictor and the outcome, without considering anyother predictors in the OLS regression model [27,28]. In this study, as thefeatures were grouped to factors, the coefficient of determination ( R ) wascalculated when each factor was in the regression, in order to quantify thecontribution of each factor.Structure coefficients quantify the direct effect of a predictor by the Pearsoncorrelation between a predictor and the predicted MPS95 value given by thefull regression model with all predictors [27,29]. In this study, the coefficientof determination ( R ) between a kinematic factor and the predicted MPS95value was calculated to quantify the contribution of each kinematic factor inthe regression.Commonality analysis (CA) [27,30] is a comprehensive method used toanalyze the relative predictive power of kinematic factors in OLS regression,as it takes into consideration not only the direct effect of each predictor, butalso the unique ( U ) and common information ( C ) when multiple predictors arepresent [30]. Commonality analysis consists of two parts: all-possible-subsetsregression and decomposition of explained variance [32]. To interpret the CAresult, the following regression sample can be considered: Y = β + β m m + β n n + (cid:15) (1)where the outcome variable Y is modeled by two predictors m and n [30]. β m , β n and β are parameters of the model and (cid:15) is regression model error. The Xianghao Zhan et al. I m pa c t s Features M PS Direction Order Power V a l ue - - . . O r de r D i r e c t i on s XY Z M A G N I T UD E P o w e r . M PS . . . . . F i g . : H e a t m a p o ff e a t u r e s i nd a t a s e t NH T S A . T h e f e a t u r e s w e r e s t a nd a r d i ze d a nd t h e v a l u e s d e n o t e h o w m u c h a f e a t u r e v a l u e d e v i a t e s f r o m t h e m e a n f e a t u r e v a l u e . T h e h e a t m a p a l s o s h o w s t h e M P S f o r e a c h s a m p l e a tt h e t o p a nd t h e d i ff e r e n t f a c t o r s o f t h e f e a t u r e o n t h e r i g h t . itle Suppressed Due to Excessive Length 7 unique information U of predictor m and n are calculated by: U m = R m,n − R n (2) U n = R m,n − R m (3)where R m,n , R n , R m are the coefficient of determination ( R ) of the model Y = β + β m m + β n n + (cid:15) and the subset models Y = β n n + (cid:15) and Y = β m m + (cid:15) respectively. The common information C of predictor m and n is calculatedby: C m,n = R n + R m − R m,n (4)The unique information of a factor can be interpreted as the R incrementsby adding the factor into a multiple regression model. The total informationof a factor is the R of the corresponding regression subset model, which canbe explained by the sum of the unique information and common information:for example, the total information of m ( T m ) is calculated by: T m = R m = U m + C m,n (5)Generally, CA decomposes the effect of each factor in a multiple regressioninto different information with R of each subset model.Dominance analysis (DA) [31] provides the analysis of dominance relation-ships among factors by calculating the increments of R when a factor is addedto every possible subset models. In DA results, level is a unique term whichdenotes the number of predictors of a subset model [27,31]. To understand theresult of DA (e.g. Table 1), level is the number of factors in the first columndenoted by k . Each entry denotes the incremental R contribution of addinga factor to a model at a specific level. Dominance of a factor in level i isdetermined by averaging all increments in R in the level i . The factor withthe largest averaged contribution dominates over other factors in this level.For example, in the DA of derivative order, the averaged contribution of thefirst-order features in level 1 ( ¯ ∆ l R o ) is calculated by:¯ ∆ l R o = 1 (cid:0) (cid:1) − R o ,o − R o ) + ( R o ,o − R o )] (6)where (cid:0) (cid:1) − R o ,o is the R ofmodel that contains zero-order features ( o
0) and first-order features ( o R o is the R of model that contains zero-order features ( o R o ,o and R o canbe interpreted in a similar manner. The averaged contribution of level 1 isrecorded in result tables where the first column is: k = 1 averaged. To avoidthe computational cost in large numbers of subset models, in the DA of powers,the powers of 1 , . , , Xianghao Zhan et al. p < .
05, see details in section 2.4), with anexception on dataset CF2 that the structured coefficients found no statisticalsignificant difference between the angular acceleration and the angular velocity( p = 0 . p < . p < .
05, see details in section 2.4). Zero-order correlation itle Suppressed Due to Excessive Length 9
MMA NASCAR NFL NHTSACF1 HMCF2 dataset pe r c en t age order MMA NASCAR NFL NHTSACF1 HMCF2 dataset pe r c en t age direction MagnitudeZYX
MMA NASCAR NFL NHTSACF1 HMCF2 dataset pe r c en t age power MMA NASCAR NFL NHTSACF1 HMCF2 dataset pe r c en t age order MMA NASCAR NFL NHTSACF1 HMCF2 dataset pe r c en t age direction MagnitudeZYX
MMA NASCAR NFL NHTSACF1 HMCF2 dataset pe r c en t age power A B CD E F
Fig. 2: Stacked bar plots of normalized mean zero-correlation coeffi-cients/structured coefficients of kinematic factors in the regression ofMPS95. The mean zero-correlation coefficients/structured coefficientswere firstly calculated in 100 iterations of bootstrapping resampling andthe results were normalized by the sum of mean zero-correlation coeffi-cients/structured coefficients of different factors. A-C. Percentage contribu-tions of features of three different derivative orders (0, 1, 2), features of fourdifferent kinematic components (x-axis, y-axis, z-axis, magnitude), featuresof seven different powers (1, 0.5, 2, 3, 4, 5, 6) given by zero-order correla-tion coefficients. D-F. Percentage contributions of features of three differentderivative orders, features of four different kinematic components, features ofseven different powers given by structure coefficients.coefficients showed the features of z-axis component were the most predictiveon dataset NHTSA ( p < . p = 0 . MM A N AS C A R C F C F H M . . . . . . . , , , , , . . . . I n c r ea s ed R Order . . . . . - .
014 0 . , , , , , . . . I n c r ea s ed R Order .
109 0 . .
067 0 . . - -
04 0 . , , , , , . . . . . I n c r ea s ed R Order . . .
015 0 . . . . , , , , , . . . . I n c r ea s ed R Order . . . . . . . , , , , , . . . . I n c r ea s ed R Order
AB E
C D F i g . : C o mm o n a li t y a n a l y s i s r e s u l t s f o r un i q u e o r c o mm o n i n f o r m a t i o n o f d i ff e r e n t p r e d i c t o r s i n t h e r e g r e ss i o n o f M P S f r o m f e a t u r e s o f t h r ee d i ff e r e n t d e r i v a t i v e o r d e r s ( , , ) o nfi v e d i ff e r e n t d a t a s e t s . ” I n c r e a s e d R ” d e n o t e s f e a t u r e s ’ c o n t r i bu - t i o n a cc o r d i n g t o c o mm o n a li t y a n a l y s i s . itle Suppressed Due to Excessive Length 11 Commonality analysis showed the predictive power varied with differentdatasets (Fig. 4) but the most predictive information always involved magni-tude. On dataset HM, the common information of features from y-axis andthe magnitude showed the highest predictive power. On dataset CF1 and CF2,the common information of features from x-axis, y-axis, z-axis and the mag-nitude was the most predictive. On dataset MMA, the common informationof features from x-axis, z-axis, and the magnitude was the most predictive.Dominance analysis showed that magnitude features dominated in the major-ity of levels of analyses with one exception on dataset MMA and NASCAR,where x-axis features dominated over other features on level 3 (when the fea-tures from y-axis, z-axis, and the magnitude were already included)(Table 1,Supplementary Tables 1-4).Table 1: Dominant analysis results for dataset HM. (a) Dominance analysis of different deriva-tive orders on dataset HM.
Order Mean R2 Additional Contribution0 1 2k=0 average k=1 average k=2 average (b) Dominance analysis of different powerson dataset HM.
Power Mean R2 Additional Contributions1 0.5 2 3k=0 average k=1 average k=2 average k=3 average (c) Dominance analysis of different components on dataset HM.
Directions Mean R2 Additional ContributionsX Y Z Magnitudek=0 average k=1 average k=2 average k=3 average .
049 0 . . . . . . - -
04 0 . . . . . . . X , Y , Z , M agn i t ude Y , Z , M agn i t ude X , Z , M agn i t ude X , Y , M agn i t ude X , Y , ZZ , M agn i t ude Y , M agn i t ude Y , Z X , M agn i t ude X , Z X , Y M agn i t ude Z YX - . . . . . . I n c r ea s ed R Direction . . .
036 0 . . . .
051 0 . .
071 0 . - . .
024 0 . . . X , Y , Z , M agn i t ude Y , Z , M agn i t ude X , Z , M agn i t ude X , Y , M agn i t ude X , Y , ZZ , M agn i t ude Y , M agn i t ude Y , Z X , M agn i t ude X , Z X , Y M agn i t ude Z YX . . . . I n c r ea s ed R Direction . . . . . . . .
005 0 . . - -
04 0 . .
178 0 .
097 0 . X , Y , Z , M agn i t ude Y , Z , M agn i t ude X , Z , M agn i t ude X , Y , M agn i t ude X , Y , ZZ , M agn i t ude Y , M agn i t ude Y , Z X , M agn i t ude X , Z X , Y M agn i t ude Z YX . . . . I n c r ea s ed R Direction MM A N AS C A R C F C F H M .
004 0 . . . . . .
016 0 . . . . .
166 0 . . . X , Y , Z , M agn i t ude Y , Z , M agn i t ude X , Z , M agn i t ude X , Y , M agn i t ude X , Y , ZZ , M agn i t ude Y , M agn i t ude Y , Z X , M agn i t ude X , Z X , Y M agn i t ude Z YX . . . . I n c r ea s ed R Direction . . . . . . . . . . - . . . . . X , Y , Z , M agn i t ude Y , Z , M agn i t ude X , Z , M agn i t ude X , Y , M agn i t ude X , Y , ZZ , M agn i t ude Y , M agn i t ude Y , Z X , M agn i t ude X , Z X , Y M agn i t ude Z YX . . . . I n c r ea s ed R Direction
AB E
C D F i g . : C o mm o n a li t y a n a l y s i s r e s u l t s f o r un i q u e o r c o mm o n i n f o r m a t i o n o f d i ff e r e n t p r e d i c t o r s i n t h e r e g r e ss i o n o f M P S f r o m f e a t u r e s o ff o u r d i ff e r e n t k i n e m a t i cc o m p o n e n t s ( x - a x i s , y - a x i s , z - a x i s a nd m ag n i t ud e ) o nfi v e d i ff e r e n t d a t a s e t s . ” I n - c r e a s e d R ” d e n o t e s f e a t u r e s ’ c o n t r i bu t i o n a cc o r d i n g t o c o mm o n a li t y a n a l y s i s . itle Suppressed Due to Excessive Length 13 p < . p > . p < . This study applied four statistical interpretation methods to analyze the re-gression of MPS95 to evaluate the predictive power of different factors of rota-tional parameters of head impact kinematics. In terms of the predictive poweramong factors when viewed individually (the direct effect defined by Budescu[31]), the zero-order correlation coefficients and structure coefficients mani-fest that the following kinematic factors generally show the highest predictivepower: 1) features based on the first derivative order of angular velocity (an-gular acceleration), 2) features based on the magnitude of the kinematics and3) features of power 1 (Fig. 2). These findings were generally supported bythe dominance analysis in which these three factors dominated in most levelsof analyses (Table 1, Supplementary Tables 1-4). The finding that magnitudefeatures were generally more predictive of MPS95 than any individual spatialdirections supports the majority of BIC, which take the magnitude of kine-matics into consideration [11,13,34]. The angular acceleration features gener-ally showed higher predictive power, which agrees with the fact that angularacceleration determines brain strain [14]. Several BIC incorporate the magni-tude of the angular acceleration into their respective computational models, . -
04 0 . . . . . -
04 6e -
04 7e -
04 0 . . . - -
04 0 . , . , , . , , , , , . , , . , , . , . , , , , . .
51 0 . . . . . I n c r ea s ed R Power MM A N AS C A R C F C F H M . . .
023 0 .
025 0 . .
004 0 . .
002 0 . . . . . . . , . , , . , , , , , . , , . , , . , . , , , , . .
51 0 . . . . I n c r ea s ed R Power . . . . . . . . .
007 0 . . . . - . . , . , , . , , , , , . , , . , , . , . , , , , . .
51 0 . . . . I n c r ea s ed R Power . . . . . . . . . . . . . . . , . , , . , , , , , . , , . , , . , . , , , , . .
51 0 . . . . I n c r ea s ed R Power . . .
009 0 . . .
004 0 . . . . . . . . . , . , , . , , , , , . , , . , , . , . , , , , . .
51 0 . . . . I n c r ea s ed R Power
AB E
C D F i g . : C o mm o n a li t y a n a l y s i s r e s u l t s f o r un i q u e o r c o mm o n i n f o r m a t i o n o f d i ff e r e n t p r e d i c t o r s i n t h e r e g r e ss i o n o f M P S f r o m f e a t u r e s o ff o u r d i ff e r e n t p o w e r s ( , . , , ) o nfi v e d i ff e r e n t d a t a s e t s . ” I n c r e a s e d R ” d e n o t e s f e a t u r e s ’ c o n t r i bu t i o n a cc o r d i n g t o c o mm o n a li t y a n a l y s i s . itle Suppressed Due to Excessive Length 15 such as RIC [13] and Power Rotation Head Injury Criterion (PRHIC) [35].The first-power features were the most predictive features on most datasets.This generally shows that the relationship between the predictor (kinematicfeatures) and the outcome (MPS95) may tend to approximate a linear rela-tionship rather than a higher order polynomial relationship.There were several dissimilarities across the datasets. First, both zero-order correlation coefficients and structure coefficients showed that angularvelocity features exhibit statistically significantly higher predictive power ondataset MMA and NASCAR (Fig. 2. The dominance analysis also supportedthe findings that the zero-order features dominated in every level on datasetMMA and NASCAR while the first-order features dominated on the otherdatasets. The higher predictive power of the angular velocity on MMA impactsmay be related to the fact that the duration of MMA impacts is generallyshorter than their counterparts in other datasets. For short-duration impacts,angular velocity tends to be more predictive of brain strain [5,14]. Based onthis finding, BrIC [12], kinematic rotational brain injury criterion (BRIC)[36], and rotational velocity change index (RVCI) [37], which all take intoconsideration the peak values of angular velocity or the numerical integral ofangular acceleration, may be more effective in predicting MPS95 on the short-duration head impacts, such as those experienced in MMA. This finding agreeswith the evaluation of BIC by Zhan et al. [15] in a recent study that foundBrIC, BRIC, RVCI, and the maximum value of angular velocity magnitudewere among the top five (out of 18 total) most accurate BIC in predictingMPS95 on MMA head impacts. On the two football datasets (CF1 and CF2),besides the most predictive magnitude features, the x-axis features and y-axisfeatures show higher predictive power than the z-axis features. This findingsuggests a potential explanation that the side or rear impacts are more likelyto lead to obvious performance decrement (OPD) of the players [26].Furthermore, on dataset NHTSA, the square-root features showed higherpredictive power than the first-power features, which were the most predic-tive on the other datasets. On dataset NASCAR, there was no statistical sig-nificance between square-root features and first-power features. This findingindicates that the underlying relationship between the features and MPS95may be slightly different among various datasets (i.e., different types of headimpacts). This may also raise concerns of interchangeably using BIC designedfor different types of head impacts, which is investigated by Zhan et al. [15].Extra model validation may be needed to show the transferability of certainBIC when applied to different types of head impacts.Additionally, commonality analysis provided more information of the pre-dictive power unique and common to different factors. In the derivative orderanalysis, the zero-order features were the most predictive factor on datasetMMA by both zero-order correlation coefficients and structure coefficients.Meanwhile, commonality analysis showed that on dataset MMA, the commoninformation of zero-order and first-order features was the most predictive,which indicates that the common information of these two factors explaineda lot of variances. These two findings collectively indicate that the subtle, unique information of the zero-order features leads the zero-order features tobe the most predictive, while the majority of predictive power is provided bythe information that are mutual to both factors. In the kinematic componentanalysis (Fig. 4), commonality analysis showed that the most predictive infor-mation was in the common information provided by y-axis features and themagnitude features for dataset HM, by x-axis, y-axis, z-axis and the magnitudefeatures for dataset CF1 and CF2, and by x-axis, z-axis and the magnitudefeatures for dataset MMA. These findings indicate the predictive informationof MPS95 may come from different spatial directions in different types of headimpacts. To accurately predict MPS95, not only the magnitude but also theother spatial directions are predictive, while the predictive power of the kine-matics component varies with types of sports. This finding indicates that thedesigns of BIC which take different spatial directions into consideration, suchas in BAM [8], BrIC [12] and RVCI [37], may bear an advantage on some typesof head impacts and provide more accurate brain injury risk evaluation.It should be noted in this study that we did not directly interpret the val-ues or the absolute values of the regression coefficients in the regression modelbecause there can be strong collinearity among the features we selected. Thecorrelated predictors in the regression model may not only lead to ambigu-ous interpretation of linearly correlated predictors but also cause a variationinflation effect [38] that leads to unstable regression coefficients. The inflatedvariance of regression coefficients can be hard to interpret as we took themean regression coefficient after bootstrapping the datasets in 100 iterationsto ensure the robustness of the analysis.Based on the findings in this study, new mathematical forms can be de-signed and validated to evaluate brain injury risk in a more accurate manner.For instance, brain injury criteria can take into consideration both angularacceleration and angular velocity [5,36], the magnitude as well as the com-ponents in different spatial directions [12,8,37], and the first-power as well asthe square-root math forms for better risk evaluation in a wider range of headimpacts.In addition to diagnostic and risk-evaluation development, upon under-standing the most predictive features in head impact kinematics, better concussion-prevention technologies can be developed with focus laid on the features pre-dictive of high brain strain. For instance, in the design process of protectiveequipment such as football helmets [39], the peaks of the angular accelerationmagnitude can be the metric to guide the design most appropriately.Although this study finds the difference in the predictive power amongfactors of head impact kinematics, there are limitations to these conclusions.Firstly, the features we analyzed are limited in power and derivative order.That is because we deemed that the higher derivative orders were highly likelyinfluenced by lower signal-to-noise ratio caused by numerical difference. Ac-cording to the present study, as the power increases when the order is higherthan 2, the predictive power generally decreases. It may be possible that thehigher derivative orders and powers show unexpected higher predictive powerin MPS95 regression. Secondly, we applied the ordinary linear regression model itle Suppressed Due to Excessive Length 17 to investigate the most direct linear relationship between predictors and theoutcome. Although several polynomial features were included to account for acertain degree of non-linearity, we did not consider other non-linearity formssuch as sinusoidal, logarithmic, exponential, etc., because any non-linear func-tion can be approximated by polynomials. It is possible that other types ofnon-linearity can better characterize the relationship between kinematics andbrain strain.This work applied four different statistical interpretation methods to ana-lyze the linear regression of MPS95 on features extracted from the kinematicsin different types of head impacts. The major findings of this work show thatthe features based on 1) angular acceleration, 2) the magnitude, and 3) thefirst power are generally shown to be the most predictive features on the ma-jority of head impact datasets. There are also dissimilarities across datasets,such as features based on angular velocity are the most predictive featureson dataset MMA and NASCAR and features based on square-root power arethe most predictive features on dataset NHTSA. The predictive informationfrom different spatial directions of kinematics varies with different types ofhead impacts. The analysis presented in this study can help the developmentof brain injury criteria and the design of impact-protection equipment. This research was supported by the Pac-12 Conference’s Student-Athlete Healthand Well-Being Initiative, the National Institutes of Health (R24NS098518),Taube Stanford Children’s Concussion Initiative and Stanford Department ofBioengineering.
References β is not enough.” Educational and Psychological Measurement61.2 (2001): 229-248.30. Ray-Mukherjee, Jayanti, et al. ”Using commonality analysis in multiple regressions: atool to decompose regression effects in the face of multicollinearity.” Methods in Ecologyand Evolution 5.4 (2014): 320-328.31. Budescu, David V. ”Dominance analysis: a new approach to the problem of relativeimportance of predictors in multiple regression.” Psychological bulletin 114.3 (1993):542.32. Prunier, J´erˆome G., et al. ”Multicollinearity in spatial genetics: separating the wheatfrom the chaff using commonality analyses.” Molecular ecology 24.2 (2015): 263-283.33. Efron, Bradley, and Robert J. Tibshirani. An introduction to the bootstrap. CRC press,1994.34. Rowson, Steven, and Stefan M. Duma. ”Brain injury prediction: assessing the com-bined probability of concussion using linear and rotational head acceleration.” Annals ofbiomedical engineering 41.5 (2013): 873-882.35. Kimpara, Hideyuki, et al. ”Head injury prediction methods based on 6 degree of freedomhead acceleration measurements during impact.” International Journal of AutomotiveEngineering 2.2 (2011): 13-19.36. Takhounts, Erik G., et al. ”Kinematic rotational brain injury criterion (BRIC).” Pro-ceedings of the 22nd enhanced safety of vehicles conference. Paper. No. 11-0263. 2011.37. Yanaoka, Toshiyuki, Yasuhiro Dokko, and Yukou Takahashi. Investigation on an injurycriterion related to traumatic brain injury primarily induced by head rotation. No. 2015-01-1439. SAE Technical Paper, 2015.38. James, Gareth, et al. An introduction to statistical learning. Vol. 112. New York:springer, 2013.39. Vahid Alizadeh, Hossein, Michael G. Fanton, August G. Domel, Gerald Grant, andDavid Camarillo. ”A Computational Study of Liquid Shock Absorption for Preventionof Traumatic Brain Injury.” Journal of Biomechanical Engineering (2020) Supplementary Table 1: Dominant analysis results for dataset CF1. (a) Dominance analysis of different deriva-tive orders on dataset CF1.
Order Mean R2 Additional Contribution0 1 2k=0 average k=1 average k=2 average (b) Dominance analysis of different powerson dataset CF1.
Power Mean R2 Additional Contributions1 0.5 2 3k=0 average k=1 average k=2 average k=3 average (c) Dominance analysis of different components on dataset CF1.
Directions Mean R2 Additional ContributionsX Y Z Magnitudek=0 average k=1 average k=2 average k=3 average a r X i v : . [ phy s i c s . b i o - ph ] F e b I m pa c t s Features M PS Direction Order Power V a l ue - - . . O r de r D i r e c t i on s XY Z M A G N I T UD E P o w e r . M PS . . . . . Supp l e m e n t a r y F i g u r e : H e a t m a p o ff e a t u r e s i nd a t a s e t H M . T h e f e a t u r e s w e r e s t a nd a r d i ze d a nd t h e v a l u e s d e n o t e h o w m u c h a f e a t u r e v a l u e d e v i a t e s f r o m t h e m e a n f e a t u r e v a l u e . T h e h e a t m a p a l s o s h o w s t h e M P S f o r e a c h s a m p l e a tt h e t o p a nd t h e d i ff e r e n t f a c t o r s o f t h e f e a t u r e o n t h e r i g h t . I m pa c t s Features M PS Direction Order Power V a l ue - - . . O r de r D i r e c t i on s XY Z M A G N I T UD E P o w e r . M PS . . . . . Supp l e m e n t a r y F i g u r e : H e a t m a p o ff e a t u r e s i nd a t a s e t C F . T h e f e a t u r e s w e r e s t a nd a r d i ze d a nd t h e v a l u e s d e n o t e h o w m u c h a f e a t u r e v a l u e d e v i a t e s f r o m t h e m e a n f e a t u r e v a l u e . T h e h e a t m a p a l s o s h o w s t h e M P S f o r e a c h s a m p l e a tt h e t o p a nd t h e d i ff e r e n t f a c t o r s o f t h e f e a t u r e o n t h e r i g h t . I m pa c t s Features M PS Direction Order Power V a l ue - - . . O r de r D i r e c t i on s XY Z M A G N I T UD E P o w e r . M PS . . . . . Supp l e m e n t a r y F i g u r e : H e a t m a p o ff e a t u r e s i nd a t a s e t MM A . T h e f e a t u r e s w e r e s t a nd a r d i ze d a nd t h e v a l u e s d e n o t e h o w m u c h a f e a t u r e v a l u e d e v i a t e s f r o m t h e m e a n f e a t u r e v a l u e . T h e h e a t m a p a l s o s h o w s t h e M P S f o r e a c h s a m p l e a tt h e t o p a nd t h e d i ff e r e n t f a c t o r s o f t h e f e a t u r e o n t h e r i g h t . I m pa c t s Features M PS Direction Order Power V a l ue - - . . O r de r D i r e c t i on s XY Z M A G N I T UD E P o w e r . M PS . . . . . Supp l e m e n t a r y F i g u r e : H e a t m a p o ff e a t u r e s i nd a t a s e t N F L . T h e f e a t u r e s w e r e s t a nd a r d i ze d a nd t h e v a l u e s d e n o t e h o w m u c h a f e a t u r e v a l u e d e v i a t e s f r o m t h e m e a n f e a t u r e v a l u e . T h e h e a t m a p a l s o s h o w s t h e M P S f o r e a c h s a m p l e a tt h e t o p a nd t h e d i ff e r e n t f a c t o r s o f t h e f e a t u r e o n t h e r i g h t . I m pa c t s Features M PS Direction Order Power V a l ue - - . . O r de r D i r e c t i on s XY Z M A G N I T UD E P o w e r . M PS . . . . . Supp l e m e n t a r y F i g u r e : H e a t m a p o ff e a t u r e s i nd a t a s e t C F . T h e f e a t u r e s w e r e s t a nd a r d i ze d a nd t h e v a l u e s d e n o t e h o w m u c h a f e a t u r e v a l u e d e v i a t e s f r o m t h e m e a n f e a t u r e v a l u e . T h e h e a t m a p a l s o s h o w s t h e M P S f o r e a c h s a m p l e a tt h e t o p a nd t h e d i ff e r e n t f a c t o r s o f t h e f e a t u r e o n t h e r i g h t . I m pa c t s Features M PS Direction Order Power V a l ue - - . . O r de r D i r e c t i on s XY Z M A G N I T UD E P o w e r . M PS . . . . . Supp l e m e n t a r y F i g u r e : H e a t m a p o ff e a t u r e s i nd a t a s e t NA S C A R . T h e f e a t u r e s w e r e s t a nd a r d i ze d a nd t h e v a l u e s d e n o t e h o w m u c h a f e a t u r e v a l u e d e v i a t e s f r o m t h e m e a n f e a t u r e v a l u e . T h e h e a t m a p a l s o s h o w s t h e M P S f o r e a c h s a m p l e a tt h e t o p a nd t h e d i ff e r e n t f a c t o r s o f t h e f e a t u r e o n t h e r i g h t . Supplementary Table 2: Dominant analysis results for dataset CF2. (a) Dominance analysis of different deriva-tive orders on dataset CF2.
Order Mean R2 Additional Contribution0 1 2k=0 average k=1 average k=2 average (b) Dominance analysis of different powerson dataset CF2.
Power Mean R2 Additional Contributions1 0.5 2 3k=0 average k=1 average k=2 average k=3 average (c) Dominance analysis of different components on dataset CF2.
Directions Mean R2 Additional ContributionsX Y Z Magnitudek=0 average k=1 average k=2 average k=3 averagek=3 average
Directions Mean R2 Additional ContributionsX Y Z Magnitudek=0 average k=1 average k=2 average k=3 averagek=3 average
Supplementary Table 3: Dominant analysis results for dataset MMA. (a) Dominance analysis of different deriva-tive orders on dataset MMA.
Order Mean R2 Additional Contribution0 1 2k=0 average k=1 average k=2 average (b) Dominance analysis of different powerson dataset MMA.
Power Mean R2 Additional Contributions1 0.5 2 3k=0 average k=1 average k=2 average k=3 average (c) Dominance analysis of different components on dataset MMA.
Directions Mean R2 Additional ContributionsX Y Z Magnitudek=0 average k=1 average k=2 average k=3 averagek=3 average
Directions Mean R2 Additional ContributionsX Y Z Magnitudek=0 average k=1 average k=2 average k=3 averagek=3 average Supplementary Table 4: Dominant analysis results for dataset NASCAR. (a) Dominance analysis of different deriva-tive orders NASCAR.
Order Mean R2 Additional Contribution0 1 2k=0 average k=1 average k=2 average (b) Dominance analysis of different powerson dataset NASCAR.
Power Mean R2 Additional Contributions1 0.5 2 3k=0 average k=1 average k=2 average k=3 average (c) Dominance analysis of different components on dataset NASCAR.
Directions Mean R2 Additional ContributionsX Y Z Magnitudek=0 average k=1 average k=2 average k=3 averagek=3 average