Brain Network Disruption in Whiplash

BACKGROUND AND PURPOSE: Whiplash-associated disorders frequently develop following motor vehicle collisions and often involve a range of cognitive and affective symptoms, though the neural correlates of the disorder are largely unknown. In this study, a sample of participants with chronic whiplash injuries were scanned by using resting-state fMRI to assess brain network changes associated with long-term outcome metrics. MATERIALS AND METHODS: Resting-state fMRI was collected for 23 participants and used to calculate network modularity, a quantitative measure of the functional segregation of brain region communities. This was analyzed for associations with whiplash-associ-ated disorder outcome metrics, including scales of neck disability, traumatic distress, depression, and pain. In addition to these clinical scales, cervical muscle fat in ﬁ ltration was quanti ﬁ ed by using Dixon fat-water imaging, which has shown promise as a bio-marker for assessing disorder severity and predicting recovery in chronic whiplash. RESULTS: An association was found between brain network structure and muscle fat in ﬁ ltration, wherein lower network modularity was associated with larger amounts of cervical muscle fat in ﬁ ltration after controlling for age, sex, body mass index, and scan motion ( t ¼(cid:2) 4.02, partial R 2 ¼ 0.49, P , .001). CONCLUSIONS: This work contributes to the existing whiplash literature by examining a sample of participants with whiplash-asso-ciated disorder by using resting-state fMRI. Less modular brain networks were found to be associated with greater amounts of cervical muscle fat in ﬁ ltration suggesting a connection between disorder severity and neurologic changes, and a potential role for neuroimaging in understanding the pathophysiology of chronic whiplash-associated disorders.

Neuroimaging findings in WAD have been mixed, with some studies reporting changes in cerebral perfusion 15,16 and white matter tract integrity, 17 and others failing to find associations by using functional or morphologic imaging. 18,19Resting-state-fMRI (rs-fMRI) is a rapidly growing tool and has been widely applied to the investigation of abnormal brain activity in clinical populations.While several studies have revealed connections between rs-fMRI and mild-to-severe traumatic brain injury, 20,21 none thus far have identified similar associations within the WAD population.
Analysis of large-scale networks by using graph theoretical approaches has recently gained traction as a method for the characterization of brain networks observed by using rs-fMRI.Within a graph theoretical framework, brain regions can be considered nodes that are linked by edges defined by the strength of pair-wise correlations between pairs of nodes.This topologic arrangement of nodes and edges can be described as a graph, which may be divided into interconnected subnetworks referred to as modules.Modularity has been found to arise in many complex systems, 22 and advances in neuroimaging have led to the characterization of brain networks as being hierarchically organized and modular systems. 23,24Modularity as a quantitative measure can be considered as the ratio between the number of connections (edges), which are located within modules, to the number of connections occurring between modules.
In this preliminary study, we investigated potential links between network modularity and symptom metrics in a sample of 23 patients with chronic WAD.6][27][28][29][30] Along with standard clinical scales for assessment of WAD, such as the Neck Disability Index 31 and the Traumatic Injury Distress Scale, 32 measures of MFI were included as an outcome metric.

Participants
This is a secondary ancillary study of participants drawn from a prospective study investigating the neuromuscular mechanisms underlying poor recovery following a whiplash injury (ClinicalTrials.govIdentifier: NCT02157038).In the main study, 97 participants were recruited (19 lost to attrition), consented, and enrolled via an urban academic emergency medicine department and were eligible provided they both reported MVC-related neck pain and were within the Quebec Task Force Classification category of WAD grade II (movement restriction with no radicular symptoms). 1Participants were eligible provided they consented and enrolled in the parent study and agreed to undergo imaging of the brain.
Volunteers were not considered if they had a spinal fracture (from the current MVC), radiologic evidence of a spinal cord lesion, or implants contraindicated for MR safety.Participants were not considered if they were pregnant or if, in the absence of an effective form of contraception, they could possibly have conceived since the first day of their last menstrual period.
The study was approved by the relevant Institutional Review Board.The 23 participants (17 female, average age of 35.8 6 12.8 years) participated in the additional data collection during their regular study visit.Participants were recruited for imaging at least 3 months following their motor vehicle collision (mean 328 6 133.16 days), and were categorized as recovered (5), moderate (10), or severe (8) in symptoms at the time of imaging.Participant characteristics are summarized in the Table.

MR Imaging Data
MR imaging data were collected by using a 3 T Prisma (Siemens) scanner with a 64-channel head-neck coil.Structural images were collected by using a 3D T1-weighted scan (TR ¼ 2.17 seconds, TE ¼ 1.69 ms, FOV ¼ 256 Â 224, 1 mm isotropic voxel size) of the brain, and a T2-weighted sagittal turbo spin-echo sequence of the cervical spine.A 3D multi-echo

Preprocessing
Muscle Fat Infiltration.MFI scores were calculated by using the Dixon water-fat scan as previously described. 12Briefly, fat and water compartments of the bilateral multifidi and semispinalis muscles from C3-C7 were manually segmented by a rater blinded to the status of the participant.The mean voxel intensity within each compartment was extracted and MFI was then calculated to generate a percentage of neck muscle fat present for each subject by using the following equation: Brain Imaging.Quality control metrics were extracted for T1 and BOLD images by using MRIQC (https://mriqc.readthedocs.io/en/stable/index.html). 33Preprocessing was accomplished by using the Nipype (https://nipype.readthedocs.io/en/latest/) 34 T1 images were bias-corrected by using ANTs (http://stnava.github.io/ANTs/)N4BiasFieldCorrection, 37 and skull stripped by using ANTs ants-BrainExtraction. Brain tissue compartments (white matter/gray matter/cerebrospinal fluid) were segmented by using FMRIB Automated Segmentation Tool (FAST; http:// fsl.fmrib.ox.ac.uk/fsl/fslwiki/fast) in native space. 38Finally, brain-extracted T1 images were transformed to the Montreal Neurological Institute 152 Nonlinear Asymmetric template version 2009c by using ANTs AntsRegistration tool. 39unctional data were section-timing corrected by using AFNI 3dTshift and motion corrected by using FMRIB Intramodal Motion Correction tool (MCFLIRT; http://fsl.fmrib.ox.ac.uk/fsl/ fslwiki/MCFLIRT) before co-registration to the T1 by using FMRIB Linear Image Registration Tool (FLIRT; https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/FLIRT). 40 The resulting transforms were concatenated with the T1-to-Montreal Neurological Institute warp and applied in 1 step by using antsApplyTransforms.Framewise displacement was used to exclude high-motion (framewise displacement .0.5) frames of the time-series, along with 1 subsequent volume.This procedure resulted in an exclusion of between 0.5%-25.8% of data, with a minimum residual time-series of 6.5 minutes (mean 9.09 6 1.5 min).To remove physiologic noise, 6 principal components of white matter/CSF signals were extracted to create aCompCor (https://nipype.readthedocs.io/en/0.13.1/ interfaces/generated/nipype.algorithms.confounds.html)nuisance regressors. 41Motion estimates along with aCompcor regressors and low frequency cosine basis regressors were removed via linear regression; regression, censoring, and smoothing (6 mm full width at half maximum) were performed by using AFNI 3dTproject.

Network Modularity Calculation
Regions of interest were defined as 5 mm-spheres by using the Power 264 coordinates. 42Temporal signal-to-noise ratios were calculated for each node/participant, and the average and standard deviation were calculated across all regions (34.0 6 10.1).ROIs were excluded if they were more than 2 SD below the mean in any participant (30 regions).Mean time course signals within the remaining ROIs were used to create a 234 Â 234 Pearson correlation matrix for each participant, which were converted to z scores.Thresholds were used to minimize the number of node pairs considered to constitute edges to a percentage of the strongest connections.To reduce the dependence of results on specific threshold selection, the connection densities used ranged from the top 2%-10% in increments of 1%, and an average across this range was taken for use in subsequent analysis, though results across specific thresholds showed similar effects (On-line Figures).This range was chosen because it was used in the construction of the Power 264 region atlas and has been utilized in previous work. 42,43he resulting undirected weighted correlation matrices were input to the Brain Connectivity Toolbox (version 2019-03-03) 44 in Matlab 2016b (MathWorks) to estimate the optimal community structure and calculate modularity.6][47] Modularity (Q) is then calculated to quantify the extent to which the network is amenable to such subdivision, with higher values of Q representing networks with a relatively high proportion of withincommunity connections to connections between communities.This form of weighted modularity was calculated as: where i and j are connections between nodes, l w is the sum of all weights in the graph, k i is the weighted degree of a node, m i is the module containing node i, and d m i ,m j ¼ 1 if m i ¼ m j , and 0 otherwise.

Statistical Analysis
Potential associations were investigated with multiple linear regression by using R version 3.4.4(http://www.r-project.org/)function lm in package stats v3.6.1. 48Covariates in all models included age, body mass index, sex (coded as 0/1 for male/female), and mean framewise-displacement (fMRI time-series motion).Clinical metrics included Neck Disability Index, 31,49 MFI, CES-D Depression Scale, 50 Traumatic Injuries Distress Scale total, Posttraumatic Diagnostic Scale hyperarousal symptom severity, 51 numeric pain rating scale (0-10), and the Hospital Anxiety and Depression Scale Depression. 52Associations were considered significant if they passed Bonferroni correction (P .05/7;P # .007).

RESULTS
Initial tests for assumptions of linear modeling did not find evidence of collinearity among predictors or extreme non-normality in clinical variables.Among models including only the covariates, the only significant effect observed was of age on MFI, with greater age being associated with larger amounts of fat infiltration (t ¼ 3.17, P ¼ .005).In the full models, network modularity was found to be negatively associated with MFI (t ¼ À4.02, partial R 2 ¼ 0.49, P , .001;Fig 1) and was not found to be associated with any other clinical metrics.Figure 2 shows the network structure for the participants with the highest and lowest MFI scores.No correlation was found between mean scan motion (framewise displacement) and modularity across participants (Pearson r ¼ À0.21, P ¼ .34).

DISCUSSION
This study is an investigation of whether rs-fMRI network measures can characterize the clinical status of a heterogeneous group of patients with persistent WAD.While previous research in WAD has yet to accurately and consistently identify markers of structural cervical spine pathology with conventional imaging, this study has found promising results of altered network structure in the brain by using more advanced imaging techniques.Such techniques have potential to influence our mechanistic understanding of WAD and other common yet enigmatic neuromusculoskeletal conditions, such as low back pain, fibromyalgia, osteoarthritis, and rotator cuff pathology.
The rationale for modularity as a chosen measure is 2-fold.First, modularity is a global metric capable of assessing wholebrain network organization without the need for defining a priori regions.This makes it an appropriate target for initial investigation given the lack of existing research into the neural correlates of WAD, as well as the modest sample size of the current ancillary study.In addition, modularity has been implicated in several  neurologic conditions that share similarities with the clinical course of WAD, primarily posttraumatic stress disorder 29 and traumatic brain injury. 53,54he prospect of anticipating the directionality of associations with modularity is not straightforward.While work in aging has suggested lower modularity may be detrimental, 55 research in psychiatric illnesses has reported deleterious effects of modularity bidirectionally, and researchers in traumatic brain injury have reported increases in the acute phase, 53 yet, the opposite in patients with persistent postconcussive syndrome. 54While some work has suggested a possible role of mild traumatic brain injury following whiplash both from kinematic modeling 56 and observations of symptom similarities, 3 the prevalence of brain injury in WAD remains largely unknown.
Given the complexities of psychosocial, traumatic, and physical components of WAD, further research is needed to determine how this finding causally relates to the condition.While there is reason to believe network changes may be related to concussive forces, it is also possible that such changes are concurrent with alterations in mood associated with trauma, 29,57 connected to symptoms of chronic pain, 58 or related to changes in a neuroimmune network, 13 in which case, the finding of disrupted modularity may not be specific to WAD but present in posttraumatic stress disorder and traumatic brain injury, for example.
The possibility of altered brain network structure in WAD raises interesting considerations for future research and clinical practice.While most patients who experience an initial whiplash injury go on to make a full recovery, many continue to exhibit symptoms for years following the event, and few quantitative tools are available for differentiating these groups in the acute stage. 59,60However, previous work has suggested cervical spine measures of MFI .20.5% (range for those with slow recovery was 6.2-40.6% and rapid recovery was 7.2-22.9%)within the first 2 weeks following the MVC resulted in a sensitivity of 87.5% (true-positive rate) and a specificity of 92.9% (true-negative rate) for predicting outcomes at 3 months post MVC. 12 The use of brain network modularity in WAD has the potential to capture a wide range of diverse network connectivity variations in a global metric that may increase prediction when used in conjunction with estimates of MFI and other clinical risk factors.
In addition to the potential diagnostic value, modularity has been shown to be predictive of treatment outcomes in contexts with potential relevance to WAD. High modularity scores have been shown to be associated with greater treatment success by using cognitive training for traumatic brain injury. 61In a similar fashion, high baseline modularity has been used to predict larger cognitive gains in response to an exercise intervention in healthy older adults. 62In light of these findings and with an ongoing trial investigating the effects of exercise on WAD (which includes fMRI), 63 network modularity presents a promising marker for treatment prediction.
This study has several limitations.As a global network measure, modularity does not carry explicit information about which sub-networks are disrupted and in what ways.This finding represents a first attempt at applying graph theoretical analysis to a sample of participants with WAD and suggests the need for additional investigation into the types of network reorganization underlying the observed change in modularity.Because no neuropsychological testing was performed, we were unable to investigate these differences in relation to cognitive outcomes such as working memory, executive function, attention, etc.
Furthermore, the intention of this work is not to suggest "we need more imaging" in clinical practice.In purest terms, judicious and informed use of advanced neuroimaging in tandem with other known risk factors may increase confidence of the primary driver of a patient's recovery trajectory, which should ultimately inform a plan of care. 64,65The work offers new directions for research in the field to consider multivariate and multisystem pre-and post-collision factors in establishing prognostic phenotypes, leading to new and more informed clinical trials.
Finally, the lack of an association between modularity and the other clinical measures is surprising and could reflect the ancillary nature of this study in that we were not powered to detect differences related to low resolution self-report measures.However, the full sample reflects the known heterogeneity of the WAD condition and provides a foundation for further mechanistic work investigating the bi-directionality of pathways linking peripheral inflammation with neural circuitries sub-serving pain, emotions, muscle structure and function, and outcomes following whiplash injury.

CONCLUSIONS
Despite evidence for the presence of cognitive symptoms, little is known about the neurobiological correlates of WAD.The discovered association between global brain network organization and a metric of WAD severity highlights the need for further advanced imaging investigations.

FIG 2 .
FIG 2. Network structure in the patients with lowest and highest MFI scores.Node colors show communities, green lines show edges within communities, red lines show edges between.The top row exhibits a high level of modularity (high within community connectivity), while the bottom row demonstrates a low level of modularity (fewer communities and more between community connections).
Note:-BMI indicates body mass index; NDI, Neck Disability Index; MFI, muscle fat infiltration; TIDS, Traumatic Injury Distress Scale; PDS, Posttraumatic Diagnostic Scale; HADS, Hospital Anxiety and Depression Scale.Â 3 mm voxel size) was collected for the quantification of MFI by acquiring images at echo times where water and fat are in phase and out of phase (TR/TE 1 /TE 2 ¼ 6.59/2.45/3.68ms) to produce water-fat ratio images.rs-fMRI was acquired with a whole-brain multiband T2*-weighted sequence (TR ¼ 613 ms, TE ¼ 22 ms, 2 mm voxel size, multiband factor ¼ 8, and 800 volumes).