Relationships between injury kinematics, neurological recovery, and pathology following concussion.

Mild traumatic brain injury affects millions of individuals annually primarily through falls, traffic collisions, or blunt trauma and can generate symptoms that persist for years. Closed-head rotational loading is the most common cause of mild traumatic brain injury and is defined by a rapid rotational acceleration of brain tissue within an intact skull. Injury kinematics-the mechanical descriptors of injury-inducing motion-explain movement of the head, which govern energy transfer, and, therefore, determine injury severity. However, the relationship between closed-head rotational injury kinematics-such as angular velocity, angular acceleration, and injury duration-and outcome after mild traumatic brain injury is not completely understood. To address this gap in knowledge, we analysed archived surgical records of 24 swine experiencing a diffuse closed-head rotational acceleration mild traumatic brain injury against 12 sham animals. Kinematics were contrasted against acute recovery outcomes, specifically apnea time, extubation time, standing time, and recovery duration. Compared to controls, animals experiencing a mild traumatic brain injury were far more likely to have apnea (P < 0.001), shorter time to extubation (P = 0.023), and longer time from extubation to standing (P = 0.006). Using least absolute shrinkage and selection operator-based regressions, kinematic parameters, including maximum negative angular velocity and time from peak angular velocity to maximum angular deceleration, were selected to explain variation in apnea time, standing time, and recovery duration. Simplified linear models employing the least absolute shrinkage and selection operator-selected variables explained a modest degree of variation in apnea time (adjusted R 2 = 0.18), standing time (adjusted R 2 = 0.19), and recovery duration (adjusted R 2 = 0.27). Neuropathology was correlated with multiple injury kinematics, with maximum angular acceleration exhibiting the strongest correlation (R 2 = 0.66). Together, these data suggest the interplay between multiple injury kinematics, including maximum negative angular velocity (immediately preceding cessation of head motion) and time from peak angular velocity to maximum angular deceleration, best explain acute recovery metrics and neuropathology after mild traumatic brain injury in swine. Future experiments that independently manipulate individual kinematic parameters could be instrumental in developing translational diagnostics for clinical mild traumatic brain injury.


Introduction
Traumatic brain injury (TBI) is a major cause of death and disability, affecting approximately 69 million individuals annually across the globe. 1 Mild TBI (mTBI) comprises the vast majority of TBIs and can induce symptoms that persist for months to years including but not limited to memory loss, sleep disruptions, emotional disturbances, headaches, fatigue, nausea, sensitivity to light, sensitivity to sound, confusion and slowed thinking. 2,3 For the vast majority of people experiencing an mTBI (frequently termed concussion), symptoms progressively abate and completely resolve within weeks or months after the injury. 4,5 However, the remaining 10-40% of people who received an mTBI will experience injury-induced symptoms that can persist for many months to years, affecting their ability to hold a job, emotional stability and quality of life. 2,4,6 Closed-head rotational injuries, the most common form of mTBI, are typically induced through rapid angular acceleration and deceleration of the head. 7-10 These movements transmit forces to the brain tissue causing diffuse damage that is difficult to detect with current clinical practices. Injury kinematics, or the mechanical descriptors of injury motion (e.g., duration of movement, angular acceleration, angular velocity, etc.), explain movement of the head and energy transfer to the brain, which determine injury severity. [11][12][13][14] However, the relationship between specific closed-head injury kinematics-such as peak angular velocity, minimum angular acceleration and positive acceleration duration-and outcomes after mTBI is currently unknown. In order to investigate this missing link, we explored the relationships between injury kinematics, acute recovery times, and neuropathology after an mTBI in swine. Within this study, we completed a retrospective analysis of archived data that were generated with a porcine model of closed-head diffuse rotational acceleration mTBI. A swine model was utilized because it can recapitulate the mechanical loading conditions observed in clinical presentations of mTBI, such as diffuse shear deformation forces which are the primary means through which most human brain injuries are generated. 10,[15][16][17] Furthermore, previous literature found that this injury model can simultaneously generate pathological and functional alterations that are observed in clinical TBI. [18][19][20][21][22][23][24][25] We postulated that acute recovery parameters including apnea time, time to extubation, return to weight bearing posture, and recovery duration would increase in animals experiencing a mTBI. Furthermore, we expected that injury kinematics would be major drivers of acute recovery parameters, and as such, that statistical models could be used to identify associations between specific parameters of injury kinematics and acute recovery.

Animal handling and anesthesia
All swine procedures were completed in accordance with the Guide for the Care and Use of Laboratory Animals 26 and followed the ARRIVE Guidelines. All protocols were approved by the University of Pennsylvania's Animal Care and Use Committee. For the current study, surgical records and injury recordings were collected from a previously completed study of adult castrated male Yucatan miniature pigs with an average weight of 34 kg. Food and water were provided ad libitum and animals were housed indoors in a facility accredited by the Association for Assessment and Accreditation of Laboratory Animal Care International.
Animals were randomly assigned to either a sham procedure (n ¼ 12) or to an injury procedure (n ¼ 24). Animals were fasted overnight prior to the injury with water remaining ad libitum. Animals were induced with a cocktail of ketamine (12-26 mg/kg) and midazolam (0.3-0.6 mg/kg). All animals were intubated with endotracheal tube and anesthesia was maintained at 1-5% isoflurane per 2-3 L of 100% O 2 for the duration of the procedure. Animals were given 0.01 mg/kg glycopyrrolate subcutaneously and eye lubricant was applied. Animals receiving an injury were given 50 mg/kg acetaminophen per rectum. Physiological monitoring of heart rate, respiratory rate and arterial oxygen saturation allowed titration of anesthesia so that all values were within acceptable ranges (heart rate between 100 and 130 beats per minute, respirations between 9 and 12 breaths per minute, and SpO 2 between 97 and 100%). A forced-air temperature management system was used to maintain normothermia throughout the procedure.
Isoflurane volume was calculated as: where G T is the total isoflurane gas provided by anesthesia equipment; G i is the percentage of isoflurane gas set by the anesthetist at each time interval; and O i is the rate of oxygen delivered through the anesthesia regulator at each time interval. The weighted anesthesia score was calculated as: where K is mass of administered ketamine; M is the mass of administered midazolam; W is the animal subject's weight; and T is the amount of time the subject was under anesthesia.
Closed-head diffuse brain injury and recovery procedure Rotational acceleration closed-head diffuse mTBI was completed while animals were under anesthesia by utilizing a HYGE pneumatic actuator. The animal's mouth was positioned around a padded bite plate and then secured to the device with adjustable snout straps. The HYGE device and accompanying linkage assembly were constructed to move the head in the coronal plane (circumferential to the brainstem) in order to produce a purely impulsive non-impact head rotation 18 (Fig. 1A). The device rapidly accelerates the swine's head and induces forces scalable to clinical TBIs. 18 Within this cohort, swine were subjected to coronal rotational injuries that ranged from 165 to 270 radians/s to induce an mTBI. Angular displacement over time was recorded with magneto-hydrodynamic sensors (Applied Technology Associates, Albuquerque, NM) connected to a National Instruments DAQ, controlled by LabVIEW. The sampling rate for the sensors was 10 kHz. Sham animals received all other procedures absent head rotation. Immediately following injury, apnea was measured as the number of seconds the animal subject did not draw breath independently. Following the injury, swine were removed from the bite plate and examined for any oral or dental injuries. Once animals were returned to their housing units, isoflurane was turned off, but oxygen continued to be delivered through the endotracheal tube. Extubation time was measured as the time from isoflurane removal to the time in which animals were extubated, as prompted by chewing on the endotracheal tube, swallowing or coughing. Standing time was measured as the time from isoflurane removal to the time in which animals were weight-bearing on all four limbs. Recovery duration was defined as the difference between standing and extubation times. Swine were continuously monitored for the duration of the recovery process. Animal data utilized in this study were from n ¼ 12 sham and n ¼ 24 mTBI archived surgical records.

Calculating injury kinematics
Angular velocity of the injury was calculated by taking the derivative of the averaged position data from the magneto-hydrodynamic sensors and was smoothed with a binomial smoothing algorithm as a low-pass filter to remove chatter. Thereafter, the derivative of the smoothed velocity trace was taken to generate angular acceleration during the injury (Fig. 1B). Six different kinematic parameters were collected from the velocity traces from each animal including: maximum angular velocity; median positive angular velocity; minimum angular velocity; time Figure 1. Schematic of injury model and resulting kinematic parameters. (A) Anesthetized swine assigned to the injury cohort were subjected to a rotational acceleration injury in the coronal plane (adapted from Ref. 18 with permission). The curved arrow indicates the rotational motion of the head during injury. (B) Movement kinematics were collected so that angular velocity and angular acceleration traces could be plotted for each injury. (C) A sample angular velocity trace and six resulting kinematic parameters describing the injury. (D) A sample acceleration trace and eight resulting kinematic parameters describing the injury. (E) Selected angular velocity and (F) angular acceleration traces from three different animal subjects exemplify differences in velocity magnitude, acceleration peaks and acceleration durations. Subject #1 (green line) is an example of a lower maximum angular velocity that occurs over a longer period of time, generating a lower, but sustained positive acceleration phase and a diminished but extended negative acceleration phase. Subject #2 (blue line) is an example of a high maximum velocity that occurs over a moderate period of time, generating a higher and sustained positive acceleration phase. Subject #3 (purple line) is an example of a high maximum angular velocity that occurs over a short period of time, generating a high but short positive acceleration phase.
to maximum velocity; time from maximum velocity to stop velocity; and time from start of injury to stop velocity (Fig. 1C). Eight different kinematic parameters were collected from the acceleration traces from each animal including: maximum angular acceleration; median positive angular acceleration; median negative angular acceleration; minimum angular acceleration; positive acceleration time; middle to minimum acceleration time; negative acceleration time; and positive and negative acceleration time (Fig. 1D).
Jerk is defined as the rate of change in angular acceleration. We were most interested in measuring the rate of change between maximum angular acceleration and minimum angular acceleration, or the maximum jerk. We calculated the maximum jerk as the largest value in the vector:J where a is angular acceleration; v is angular velocity; p is position; and t is time. Excursion was defined as the distance traveled during the injury in radians. We calculated excursion at maximum angular velocity as: where p V¼max is position at maximum angular velocity and p V¼start is position at the start of the injury. We calculated total excursion of the injury as: where p V¼stop is position at the end of the injury. All kinematic traces were processed and individual parameters were calculated in MATLAB version 2019b (9.7.0.1247435).

Sacrifice and tissue acquisition
At the designated time point, animals were induced and intubated as described above. Thereafter, animals were transcardially perfused with 0.9% heparinized saline followed by 10% neutral buffered formalin. Decapitated tissue was stored overnight in 10% neutral buffered formalin before the brain was extracted and post-fixed in 10% neutral buffered formalin for 1 week at 4 C. Brains were blocked coronally every 5 mm, paraffin embedded and 8 mm sections were collected via rotary microtome.

Pathological characterization
Assessment of diffuse axonal injury, a key characteristic of TBI pathology, was previously quantified and reported for this animal cohort. 27 Briefly, the extent of amyloid precursor protein (APP) was measured by an analyst blinded to injury condition on an ordinal scale from 0 to 3 with 0 representing no pathology and 3 representing severe pathology. Because the extent of APP pathology changes over time after the injury, we only considered pathological scores from animals sacrificed 7 days after injury. Pathology scores were available for periventricular white matter, striatum, ventral thalamus, dorsal thalamus, hippocampus/fornix and cerebellum. The total APP pathology score was generated by taking the mean pathology score across these six regions. The correlative relationships between pathology and several injury kinematics were described with a line of best fit and 95% confidence interval bands.

Statistical analyses
As described above, archived records from swine surgical procedures yielded information related to the dose and duration of anesthesia, as well as the time to recovery after a sham or mTBI procedure. Twenty-four injured and 12 sham animals meeting the inclusion criteria of having complete and interpretable surgical records were utilized. As none of the sham animals exhibited apnea, we categorized apnea as present or absent and used Fisher's exact test to assess differences in the proportion of animals with apnea. To compare sham and injured animals, Welch two sample t-tests were utilized to assess mean differences in extubation time, standing time, and recovery duration. Pairwise correlations were visualized to investigate associations between outcome measures and anesthesia for all animals, between outcome measures and kinematics for all animals, and between outcome measures and kinematics for mTBI animals.
Principal component analysis (PCA) was used for dimensionality reduction within the mTBI animal cohort. To stabilize the model, we used only kinematic variables with a pairwise Spearman's correlation coefficient jqj < 0.9, leaving 12 kinematics as input variables for PCA, after Z-scoring. Importantly, PCA was completed with the kinematic variables, but not the recovery measures. The principle components were plotted in a bivariate fashion and color-coded to visualize which components, and thus which combination of kinematics, were most strongly associated with recovery outcomes. Then, linear models were fit for each outcome metric with the first four principal components (PCs) as additive input variables. Because PC2 (the second PC) explained the largest proportion of variance, we generated scatter plots of PC2 values against apnea time, extubation time, standing time, and recovery duration.
Next, we used least absolute shrinkage and selection operator (LASSO) regressions to select individual kinematic variables associated with outcomes. 28 The same 12 kinematic variables were used as candidates. LASSO is a form of penalized regression (see Supplementary methods). Given a penalty and a set of candidate predictor variables, the method optimizes a function such that the coefficients of certain variables are set to zero. The penalty was determined using 10-fold cross-validation. Specifically, a candidate set of 100 penalties was pre-specified, and for each penalty the model was fit to the data for 22 randomly chosen subjects. The fitted models were used to predict outcomes for the remaining 2 subjects and the mean square error determined for the prediction set (2 subjects). This process was repeated 10 times; the cross-validated penalty minimized the average mean square error across the 10 prediction sets. Then, using the complete dataset (24 subjects) we fit the LASSO regression model using the cross-validated penalty and reported variables with nonzero coefficients. Lastly, to assess the potential of the variables to explain variation in outcome, we refit the model using only those variables selected by the LASSO regression, and reported adjusted coefficients of variation R 2 . Because these original data were also used for variable selection, these adjusted R 2 values are expected to overestimate the predictive value of the model in a new dataset.
We further investigated the stability of our LASSOregression results, using leave-one-out paradigm. Specifically, the LASSO was fit to a subset of 23 animals (training set) and the variables with non-zero coefficients were recorded. This process was repeated 24 times so that each animal was left out once. The likelihood of each kinematic variable remaining in the LASSO model was calculated as a percentage of these 24 iterations. The models from the training sets were applied to the left-out animal (prediction set) and the predicted values from the training model were plotted against the observed values.
For all tests, the type I error rate was set to 0.05 and all hypothesis tests were two-sided. No adjustments for multiple comparisons were performed for the individual predictors. All statistical analyses were carried out in R Studio Version 3.6.3 using R and the Hmisc, stats, rstatix, corrplot, Matrix, glmnet, factoextra and ggplot2 packages. Line graphs, column graphs and scatter plots were generated in Prism version 8.4.3(471).

Data availability
The datasets used for the current study are available from the corresponding authors upon reasonable request.

Closed-head diffuse TBI influences recovery outcomes
Mild TBI was generated in male Yucatan miniature swine with an inertial closed-head diffuse brain injury in the coronal plane. Injuries were classified as mild TBIs because evidence of gross pathology, overt bleeds or structural changes in MRI (data not shown) or following tissue explantation were extremely rare (e.g., images of gross structure and histology please see Grovola et al. 29 ). Angular velocity and angular acceleration were plotted over time for each injury and 17 resulting kinematic parameters were collected to describe the injury (Fig. 1).
Deployment of the HYGE injury device allows for control over the maximum angular velocity experienced by the animal subject, defined as the peak value in the angular velocity trace (Fig. 1C). As a result, injured animals in this cohort intentionally experienced coronal TBIs ranging from 165 to 270 rad/s. However, maximum angular velocity is not modulated in isolation. In general, when maximum angular velocity decreases, maximum angular acceleration decreases, minimum angular acceleration (e.g. maximum deceleration) increases, the duration of the positive acceleration phase increases, and the duration of the negative acceleration phase increases ( Fig. 1E and F). After mTBI, none of the control animals experienced apnea, while 80% of the injured animals experienced apnea that ranged from 7 to 44 s (P < 0.001; 95% CI ¼ 0.000-0.146). Contrary to our expectations, mean time to extubation declined by 4.58 min (P ¼ 0.023; 95% CI ¼ 0.706-8.461) in injured animals relative to sham animals. The mean time to standing was extended by 12.71 min in injured animals although this trend was not significant (P ¼ 0.063; 95% CI ¼ À0.753 to 26.169). The mean recovery duration was extended by 17.29 min in injured relative to non-injured animals (P ¼ 0.006; 95% CI ¼ 5.382-29.202) (Fig. 2).

Injury kinematics are correlated with outcome measures
We next explored how different injury parameters related to one another using pairwise Spearman correlations. In the visual display, the upper half of the matrix depicts the correlative relationship with an ellipse while the bottom half shows the actual estimate. The shape and color of the ellipse in each cell denote the magnitude and direction of the Spearman's correlational coefficient (q) for the factor listed above and to the left. A correlation matrix with the four outcome parameters (black text) and the five anesthesia parameters (grey text) utilizing data from sham and injured animals appears in Supplementary Fig. 1. We did not observe strong correlations between anesthesia parameters and recovery outcomes. Ketamine dose was modestly correlated with apnea time, standing time and recovery duration (jqj 0.36). The weighted anesthesia score, which factors in normalized ketamine, normalized midazolam, anesthesia time and isoflurane volume was weakly correlated with each of the four outcome parameters ( Supplementary Fig. 1).
Next, we generated a correlational matrix between kinematic parameters and recovery outcomes for injured animals (Fig. 3). These data suggest the pairwise correlation between outcome measures and injury kinematics was much stronger than the correlation between outcome measures and anesthesia parameters. Among the outcome variables, apnea time was positively correlated with extubation time but was largely independent of the standing time and recovery duration. Extubation time was weakly correlated with standing time but was not correlated with The correlation matrix depicts the Spearman's correlation coefficient (q) for each combination of injury parameters. Only animals experiencing an mTBI were utilized in this matrix, no sham animals were included. The matrix is organized so that the upper half of the matrix depicts the magnitude and direction of Spearman's correlational coefficient with coloured ellipses while the lower half of the matrix listed the calculated Spearman's correlational coefficient. The orientation and color of each ellipse is indicative of the magnitude and direction of the correlation with navy ovals indicating a strong positive correlation, pastel or white circles indicating weak or no correlation, and burgundy ovals indicating a strong negative correlation. Injury parameters were colored to show outcome measures (black text), kinematic parameters with positive values (grey text), or kinematic parameters with negative values (red text). All references to 'velocity' and 'acceleration' refer to 'angular velocity' and 'angular acceleration', respectively. (D) The time from extubation to the time in which animals were weight-bearing on all four limbs was significantly increased between sham and mTBI animals (P ¼ 0.006). Data represent mean 6 standard deviation with *P < 0.05; **P < 0.01; and ***P < 0.0001. recovery duration. Standing time and recovery duration were highly correlated due to their mathematical relationship (Fig. 3). We observed that apnea was only moderately correlated with minimum angular velocity and minimum angular acceleration (jqj ! 0.40) and both of these relationships were negatively correlated. Similarly, extubation time was only correlated with minimum angular velocity (jqj ! 0.40) and did not exhibit strong relationships with any other kinematic terms. In line with our finding that mTBI reduced extubation time, we observed negative correlations between extubation time and the majority of the injury kinematics. Standing time and recovery duration were correlated with many kinematic terms. Standing time had jqj ! 0.40 for six kinematic variables while recovery duration had jqj ! 0.40 for seven kinematic variables. Because recovery duration includes standing time in its calculation, it was not surprising to find that standing time and recovery duration exhibited very similar relationships to many of the kinematic terms (Fig. 3). Time to maximum velocity, median positive angular acceleration and positive acceleration time were the three most strongly correlated kinematics for standing time (jqj ! 0.50). Median positive angular acceleration and positive acceleration time were the most strongly correlated kinematics for recovery duration (jqj ! 0.50).
We generated a similar correlation matrix between outcome metrics and kinematic parameters for all sham and injured animals ( Supplementary Fig. 2). Because none of the sham animals exhibited apnea, there was a point mass at the origin which artificially inflated the relationships between apnea and kinematics. Several scatter plots represented within this figure are individually plotted in Supplementary Fig. 3.

PCA separates outcomes according to injury kinematics
Outcome metrics following mTBI were highly variable. Specifically, variability in apnea time and recovery duration suggest that factors other than just the presence or absence of injury may contribute to the spread of the data. As previously mentioned, animals within the mTBI cohort experienced vastly different injuries from one another. Indeed, maximum angular velocity ranged from 165 to 270 rad/s, maximum angular acceleration ranged from 36 004 to 126 237 rad/s 2 , and minimum angular acceleration ranged from À66 600 to À185 700 rad/s 2 ( Supplementary Fig. 4). The pairwise correlations provide information about univariate associations between kinematic variables and outcome. However, we wanted to know if some combination of kinematic parameters could explain acute recovery outcomes. Therefore, we completed PCA in order to reduce the dimensionality of the mTBI animal kinematic dataset. PCA was completed using 12 input kinematics that were not collinear with each other and PCA did not include outcome parameters in the analysis to ensure that any discovered relationships would not be exaggerated due to overfitting the data.
The first four PCs had eigenvalues greater than 1, indicating that they each explained more variance than any single kinematic input parameter (Fig. 4A). Together, the first four PCs explained 93.7% of the variance (Fig. 4B). PC1 explained 45.9%, PC2 explained 27.6%, PC3 explained 11.6% and PC4 explained 8.6% of the variance. A variable map of PC1 versus PC2 illustrates the distribution and weighting of the variables across the different PCs (Fig. 4C).
We plotted PC1 against PC2 and colored the points according to apnea time, extubation time, standing time or recovery duration (Fig. 4D-G). Segregation between high and low recovery outcomes was not obvious but modest separation did occur. Animals experiencing shorter apnea times generally clustered in lower values along PC2 (Fig. 4D). PCA did not separate animals according to extubation times (Fig. 4E). Animals experiencing shorter standing times clustered in PC2's higher values and PC1's lower values while animals experiencing longer standing times clustered in PC2's lower values (Fig. 4F). Animals experiencing shorter recovery durations clustered in PC2's higher values and PC1's lower values while animals experiencing longer recovery durations clustered in PC2's lower values (Fig. 4G).
Next, we generated linear models of the sum of the first four PCs in order to determine which PCs were most strongly related to each outcome metric (outcome $ PC1 þ PC2 þ PC3 þ PC4). We determined that PC2 was the only PC that was significantly associated with apnea times (P ¼ 0.042), standing time (P ¼ 0.014), and recovery duration (P ¼ 0.004). None of the PCs were significantly associated with extubation time. Because PC2 values were most strongly associated with the majority of the outcomes, we plotted the PC2 value for each animal against apnea time, extubation time, standing time, and recovery duration (Fig. 5). PC2 was weakly related to apnea time (R 2 ¼ 0. 19), was not at all related to extubation time (R 2 ¼ 0.03), and was modestly correlated to standing time (R 2 ¼ 0.26) and recovery duration (R 2 ¼ 0.34) (Fig. 5). The variables that most strongly contributed to PC2 were minimum angular velocity, excursion at maximum velocity, minimum angular acceleration, negative acceleration time, total excursion, and positive acceleration time.

Identified injury kinematics jointly explain acute recovery after TBI
Our PCA results suggest that some combination of injury kinematics can partially explain acute recovery outcomes. Building on these findings, we used LASSO regression to screen for and select kinematic variables relevant to acute recovery outcomes. Like the PCA, we utilized the same 12 kinematic terms that were not collinear. In using the LASSO regression, we determined that the variation in apnea time was partially explained by a linear model of three terms: minimum angular velocity, minimum angular acceleration, and middle to minimum acceleration time (Fig. 6A). Variation in standing time and recovery duration were explained by linear models with six terms: maximum angular velocity, minimum angular velocity, excursion at maximum velocity, positive acceleration time, negative acceleration time, and middle to minimum acceleration time (Fig. 6A). Two kinematic terms-minimum angular velocity and middle to minimum acceleration time-were preserved across the apnea time, standing time, and recovery duration models. The LASSO regression eliminated all the input kinematic variables when generating the model for extubation time, suggesting that no variables had sufficiently strong associations to remain in the model. After LASSO selected kinematic terms for each recovery parameter, we refit the data with a linear model in order to determine the best linear fit using only the LASSOselected coefficients. We observed that these linear models could modestly relate kinematics to apnea (adjusted R 2 ¼ 0.18), standing time (adjusted R 2 ¼ 0. 19), and recovery duration (adjusted R 2 ¼ 0.27). To understand the accuracy of these three equations, we next calculated the apnea time, standing time, and recovery duration for all injured animals by plugging in each animal's kinematic values into each equation. These calculated values were plotted against the actual, experimental outcome values. The LASSO calculated values for apnea exhibited a modest correlation with the actual apnea time (Fig. 6B). The LASSO calculated values for standing time exhibited a modest positive correlation with the actual standing time (Fig. 6C). The LASSO calculated values for recovery duration also exhibited a moderate positive correlation with the actual recovery duration (Fig. 6D).
The LASSO regression was repeated on z-scored kinematic data so that the coefficients could be compared using the same scale. For apnea, changes in minimum angular velocity were associated with the largest changes in outcome followed by middle to minimum acceleration time (Supplementary Table 1). For standing time, changes in negative acceleration time were associated with the largest change in outcomes. For recovery duration, changes in maximum angular velocity were associated with the largest change in outcome (Supplementary Table 1).
While our study is unique in terms of analyzing animal behavior following TBI, with a sample of only 24 animals, it seemed unlikely that we could effectively create and assess a predictive model by splitting the data into a training and test set. 30 However, to assess the robustness of our variable selection process, we calculated the probability that each kinematic term would be preserved in each recovery equation by iteratively generating LASSO models on a smaller training cohort. The models for apnea time frequently included minimum angular velocity (92%), minimum angular acceleration (79%), and middle to minimum acceleration time (83%) (Supplementary Fig.  5A). The recovery duration equation had several terms that were consistently expressed in the model including maximum angular velocity (88%), minimum angular velocity (88%), and middle to minimum acceleration time (79%). The standing time equation exhibited coefficients with much lower likelihoods of incorporation such as excursion at maximum velocity (54%) (Supplementary Fig.  5A). The terms maximum angular acceleration, maximum jerk, median negative angular acceleration, and positive and negative acceleration time were rarely or never used within any of the recovery models ( Supplementary Fig.  5A). Not surprisingly, given the sample size, the leaveone-out samples were poorly predicted by the models using only 23 animals ( Supplementary Fig. 5B-D).
Brain pathology is correlated to injury kinematic parameters APP accumulation in axons is the gold-standard used to identify diffuse axonal pathology following closed-head TBI. [31][32][33][34] Previous efforts in our lab have characterized the distribution of APP brain pathology in these animals over time. 27 However, these pathology metrics have not been compared against injury kinematics. We collected pathology data from animals surviving seven days after a sham or coronal injury. Average APP pathology was defined as the mean pathology score of 0 (no pathology) to 3 (severe pathology) across six brain regions: periventricular white matter, striatum, ventral thalamus, dorsal thalamus, fimbria/fornix, and cerebellum. We plotted the average APP pathology score against eight kinematic parameters that were important in describing the injury. We plotted average APP pathology against terms that described the magnitude of the injury: maximum angular velocity, maximum angular acceleration, minimum angular acceleration, and maximum jerk. We also plotted average APP against terms that were deemed important in the LASSO and PCA analyses: minimum angular velocity, middle to minimum acceleration time, median negative angular acceleration, and excursion at maximum velocity. We observed that the average APP pathology score was correlated (R 2 > 0.45) with maximum angular velocity, maximum angular acceleration, minimum angular acceleration, maximum jerk, middle to minimum acceleration time, and median negative angular acceleration (Fig. 7A-F). The average APP pathology score was moderately correlated with excursion at maximum velocity and had no correlation with minimum angular velocity (Fig. 7G-H).
Maximum angular acceleration was the most strongly correlated with average APP pathology scores with R 2 ¼ 0.66.
We also investigated how average APP pathology correlated to our four recovery parameters: apnea time, extubation time, standing time, and recovery duration. Correlations between pathology at 7 days and apnea time, extubation time, and standing time were weakly associated, as demonstrated by R 2 values <0.30 in all cases ( Supplementary Fig. 6A-C). Pathology at 7 days and recovery duration was moderately correlated with an R 2 value of 0.40 ( Supplementary Fig. 6D).

Discussion
Diffuse closed-head rotational acceleration is the most common form of TBI, 7-10 but the specific injury kinematic parameters of head rotation that drive pathology and affect neurological recovery are currently unknown. This retrospective study shows that acute recovery parameters including apnea time, extubation time, and recovery duration are altered after a single mTBI in swine and that a number of kinematic parameters are associated with recovery. We also found that injury kinematic parameters of closed-head diffuse injury are correlated with neuropathological outcomes after mTBI. To our knowledge, this study represents the first attempt to understand the contributions of distinct injury kinematics to recovery using a large animal model of TBI. Taken together, these data inform how injury kinematic parameters affect acute outcomes and suggest that understanding the interplay between multiple injury kinematics may be key in predicting recovery following trauma. A major goal of this work was to determine if combinations of injury kinematics could explain variation in recovery outcomes. The first step in investigating the multicomponent relationships was to complete PCA. We observed that of all the PCs, PC2 was most strongly associated with the recovery outcomes apnea, standing time, and recovery duration. PC2 factored in all 12 of the input kinematics but minimum angular velocity, excursion at maximum velocity, minimum angular acceleration, negative acceleration time, total excursion, and positive acceleration time were the factors that most strongly affected PC2. PCA data suggest that recovery outcomes could be better explained with a combination of kinematics instead of any singular term.
To build on this PCA finding and to improve the interpretability of the statistical relationships to recovery parameters, we employed a LASSO regression, a variable screening and selection method. The results of the LASSO suggested a modest component of the variance in outcome was explained by the kinematic parameters and identified specific parameters of interest for further mechanistic studies. This regression analysis determined that minimum angular velocity, minimum angular acceleration, and middle to minimum acceleration time were the most important factors in explaining apnea time. While all three kinematics contribute to the model's accuracy, minimum angular velocity was the most strongly weighted factor. This sparse and stable linear model is consistent with the observation that apnea time was only strongly correlated with a few kinematic terms in the correlation matrix. In contrast to the apnea equation, maximum angular velocity, minimum angular velocity, excursion at maximum velocity, positive acceleration time, negative acceleration time, and middle to minimum acceleration time were the most important factors in explaining both standing time and recovery duration. These models preserved the same six kinematic inputs, which is probably due to the mathematical relationship that standing time and recovery duration share. Within the standing time and recovery duration equations, negative acceleration time and maximum angular velocity were the most strongly weighted factors, respectively. Minimum angular velocity and middle to minimum acceleration time were the only terms that were preserved in all three equations suggesting that these kinematic variables play a leading role in determining several acute recovery outcomes.
We were surprised to find that minimum angular velocity (i.e., maximum negative angular velocity) played a leading role in the PCA analysis and was selected by the LASSO regression in each model. Indeed, minimum angular velocity is a kinematic term that does not have a large magnitude and is not associated with the largest changes in angular acceleration. The minimum angular velocity term represents a brief change in direction of the HYGE's linkage arm's movement at the end of the excursion, immediately prior to cessation of head motion. While the HYGE device is pneumatically driven, the acceleration and deceleration characteristics of the actuating piston are determined by the geometry of acceleration and decelerating metering pins (i.e., cones surrounding the piston shaft), respectively. In particular, these internal components determine the characteristics of deceleration-and eventual stopping-of the piston, rather than the piston having a 'hard stop' which would occur if motion ceased due to the piston or the side arms striking a rigid surface. As such, the large negative angular acceleration (i.e., rapid deceleration) preceding the stoppage of motion of the HYGE's linkage arm typically generates a small bounce of the animals' head when it reaches the stopping position. That this small change in directionand perhaps contemporaneously occurring reversals in angular acceleration and jerk-could be one of the larger drivers of acute recovery outcomes, may speak to the complexity of viscoelastic tissue strains occurring during rapid transitions between rotational acceleration and deceleration. Indeed, because we measured kinematics from devices fixed to the HYGE's linkage arm, not the brain directly, small changes in minimum angular velocity for the HYGE's linkage arm may insinuate larger changes in movement inside the intracranial space. Future work attempting to control the magnitude and duration of the minimum angular velocity 'bounce' of the HYGE's linkage arm-as well as resulting effects on head and brain tissue motion-would help better understand the contribution of this kinematic parameter to outcomes.
We also reported correlational relationships between injury kinematics and white matter pathology 7 days after the injury. Unlike the PCA and LASSO analyses, maximum angular acceleration, maximum angular velocity, and median negative angular acceleration had the strongest relationships with pathological outcome. Minimum angular velocity, which played a leading role in recovery terms, had no correlation to pathology. Middle to minimum acceleration time, another lead kinematic in explaining recovery, had a moderate correlation to pathology, although it was weaker than several other kinematics. This may suggest that the kinematics that drive acute recovery parameters are not necessarily the same kinematics that contribute to white matter pathology. This also suggests that studies utilizing closed-head rotational models should report angular velocity kinematics, angular acceleration kinematics, and injury timing information in order for others to replicate their work. Future studies investigating why white matter pathology at 7 days did not directly relate to acute recovery parameters will also be helpful in contextualizing these findings.
Implications for developing explanatory relationships between injury kinematics and recovery outcomes are immense. Injury kinematic data can be utilized to develop computational models of brain movement during injury. [35][36][37] Many kinematic parameters are related to one another and computational modelling may be able to enhance predictive accuracy by informing which parameters to prioritize. If head kinematics could be measured as the injury occurs (through mounted sensors in civilian cars or on military/athlete helmets), then physicians may be able to utilize the kinematic parameters in a computational model to understand the injury severity and predict the patient's recovery experience. 14,[38][39][40][41][42] Indeed, predictive models could be especially useful for mTBI because no overt bleeding or pathology may be evident even though many patients report persistent effects for months or years after an injury. [2][3][4] In line with these goals, we completed injuries with kinematics that mimic clinical presentations of mTBI or concussion. Clinical concussion is associated with confusion, headache or dizziness; loss of consciousness that does not exceed 30 min; loss of memory that lasts <24 h; and no gross structural changes. 43 We did not have any evidence of brain bleeds in these animals in MRI scans (data not shown) or in the tissue during explantation and subsequent histopathological analyses. 29 While it is difficult to extrapolate clinical characteristics of consciousness and amnesia to subjects that cannot self-report, 44 we observed that recovery duration was extended by 17 min on average in injured swine relative to sham animals. It took approximately 25 min for sham animals to recover from the anesthetic effects after inhalants were removed. Extending recovery duration by 17 min fits into the 30 min threshold that is used for clinical definitions of concussion. In addition to utilizing gross pathology and recovery physiology to determine injury severity, injury forces experienced by swine can be scaled to clinically relevant forces experienced in human brain during TBI by utilizing Holburn's scaling equation. 18 Within this study, swine experienced minimum angular acceleration values that ranged from À66 000 to À186 000 rad/s 2 with an average minimum angular acceleration of À128 000 rad/s 2 . Scaling these according to the brain masses of swine and humans, we estimate that this would be approximately equivalent to peak angular acceleration values of À10 000 to À29 000 rad/s 2 in a human brain. These injury levels are proportional to some collision events experienced by high school and collegiate football players ( 15 000 rad/s 2 ) 13,45 or motor vehicle crash simulations ($25 000 rad/s 2 ). 46 However, it is important to note that comparing pig rotational acceleration values to clinical events is challenging because real-world scenarios will not conserve the distance between the center of mass and the center of rotation, which would affect the amount of force experienced. Furthermore, clinical cases of TBI will not confine movement to one rotational plane, as we have in this experiment.
This study was based on a retrospective analysis of archived injury kinematics and surgical records from swine experiencing either a sham or an mTBI. To our knowledge, this study represents one of the largest studies of the behavioral effects of TBI in swine. With that said, a study of a complex behavioral outcome in only 24 subjects, was likely underpowered to detect more subtle associations. Given the sample size, splitting the sample into a training and test sample for the purpose of predictive modelling was not warranted. Future behavioral studies are needed to determine the reproducibility of our findings. In addition, as is the nature of retrospective studies, we only utilized outcome measures that were consistently recorded across all animals. However, future studies correlating injury kinematics to long-term behavioral and physiological outcomes could be particularly informative to clinical translation of these models. Indeed, we suspect that independently modulating injury kinematics could have drastically different outcomes on pathology, physiology, and behavior. 10,36 Because the brain is viscoelastic, the duration over which kinematics are applied to the brain determines how deep injury forces penetrate into the tissue, with rapid injuries concentrating forces on the surface of the brain and with extended injuries distributing forces into deep structures of the brain. Therefore, we suspect that longer injury durations may generate more neuronal damage, synapse disruption, and gliosis in deeper brain structures, potentially affecting consciousness, mood, addiction and attention. 10,[47][48][49] In contrast, we speculate that injuries with higher maximum angular velocity and maximum jerk may generate more neuronal damage, synapse disruption, and gliosis in cortical brain structures, potentially affecting aggression, memory and coordination. 48,50 These results illustrate that there is a need for future studies to independently control each kinematic parameter during TBI to fully understand their contributions to pathology, physiology and behavior.
Additionally, as inertial models of closed-head diffuse TBI are becoming more common in the neurotrauma field, there is a need to compare and contrast results across institutions and investigators. 51,52 Indeed, several research groups employ a similar porcine injury model of diffuse closed-head mTBI because the architecture of swine's large gyrencephalic brains allows loading conditions that are scalable to clinical TBIs. 21,51-54 Our findings suggest that maximum angular velocity alone may not drive recovery parameters or neuropathology. Systematically controlling and reporting kinematic parameters that influence pathology and outcome will be critical to compare and contrast findings between studies.
In completing this study, we expected that an injury would delay all recovery outcome parameters. However, in contrast to our expectations, we found extubation times were significantly decreased in injured animals relative to sham animals. Future studies investigating the causal relationship between injury and extubation time are necessary to tease apart this unexpected relationship.
Moreover, coupling recovery outcomes of apnea time, standing time, and recovery duration with physiological changes in the brain could be informative in interpreting these data. Indeed, investigating how pain, stress and brainstem pathology affect extubation and standing times could support this work. Building on this, apnea time is not commonly reported in closed-head rotation models of swine TBI. Investigating potential causes of apnea and the implications of apnea (transient reduction in oxygenation) for other behaviors after injury are important in contextualizing these findings. Furthermore, contrary to our expectations, we found that anesthesia factors-within the ranges administered in this study-did not play a major role in the outcome parameters.
While the results of this study suggest that kinematic parameters can be used in combination to relate apnea time, standing time, and recovery duration; this study has several notable limitations. Most importantly, this retrospective study did not allow independent modulation of individual injury kinematics. Here, statistical analyses attempt to weight individual kinematic parameters in order to explain outcomes. However, independently controlling one kinematic parameter at a time will facilitate a better understanding into how that type of injury affects pathology, physiology, and animal behavior. Additionally, explanatory models in this study utilized linear models with additive relationships. Modelling outcomes with non-linear relationships or investigating interactions of kinematics could enhance the accuracy of these models but requires much larger sample sizes for statistical validity. In addition to these limitations, the retrospective nature of this study underscores the importance of consistent record keeping related to post-injury acute recovery metrics, thus mitigating challenges in interpreting archived procedure records. Animal recovery progress was described approximately every 10 min until an hour after swine were completely ambulatory. Only complete records that had descriptive narratives of acute recovery were included in the study. Because these acute recovery parameters were only deemed a dependent variable after all surgeries had been completed, we suspect that inaccuracies and variability would be conserved across sham and injured animals. Furthermore, personnel involved in the study would have been unaware that acute recovery times would be a part of this study and so would have acted in an unbiased manner. In this way, the retrospective study was limited in temporal resolution and accuracy, but also had the benefit of decreased observer bias. Future studies in which kinematic variables are independently modulated and dependent variables are systematically recorded could be informative in validating these results.
In addition to limitations of the study design, we were only able to correlate injury kinematics with acute metrics of recovery. However, many other parameters could be systematically quantified in swine after injury. Future swine TBI studies could record other acute outcomes including extent of vocalization, control of head movements, supporting weight on two limbs, return of coordinated walking, gait rigidity, rearing, eye movements, and oculocephalic reflexes. 55 Furthermore, the development of long-term behavioral metrics for swine to test learning, memory, aggression, anxiety, saccades, sleep disturbances, and addiction could act as a critical link in understanding injury side-effects in clinical presentations. In addition to investigating how these injuries affect behavioral outcomes, altering the plane of rotation would be critical for translation to clinical predictive models. Typically, clinical TBIs are not constrained to closed-head rotational acceleration within one plane, but rather movement of the head occurs in multiple planes during an injury. Indeed, previous work in pediatric swine suggests that TBIs generated in the sagittal plane may generate more severe neuropathology, increase the duration of unconsciousness, increase the likelihood of apnea, and reduced animal subjects' activity levels compared with TBIs generated in other planes. 56,57 Investigating how kinematics affect recovery in injuries generated in the sagittal plane, axial plane, or oblique rotations would be important to combine with this work to predict clinical recovery. Finally, sex was not accounted for in this study, but it would be interesting to see if these trends were preserved across both male and female swine.
In summary, we utilized a large-animal model of closed-head rotational TBI that can scale loading conditions to TBI in humans. We found that acute recovery after mild closed-head diffuse TBI in swine could be modestly described with an additive linear model informed by injury kinematics. Together, these data suggest that no singular kinematic parameter drove outcome, rather, interplay between multiple injury kinematics contribute to recovery parameters and neuropathology after mTBI in swine. We believe that better understanding of how injury kinematics affect neurological recovery and neuropathology will inform the development of more translational diagnostic criteria for TBI in the future.