The "vestibular neuromatrix": A proposed, expanded vestibular network from graph theory in post-concussive vestibular dysfunction.

Convergent clinical and neuroimaging evidence suggests that higher vestibular function is subserved by a distributed network including visuospatial, cognitive-affective, proprioceptive, and integrative brain regions. Clinical vestibular syndromes may perturb this network, resulting in deficits across a variety of functional domains. Here, we leverage structural and functional neuroimaging to characterize this extended network in healthy control participants and patients with post-concussive vestibular dysfunction (PCVD). Then, 27 healthy control subjects (15 females) and 18 patients with subacute PCVD (12 female) were selected for participation. Eighty-two regions of interest (network nodes) were identified based on previous publications, group-wise differences in BOLD signal amplitude and connectivity, and multivariate pattern analysis on affective tests. Group-specific "core" networks, as well as a "consensus" network comprised of connections common to all participants, were then generated based on probabilistic tractography and functional connectivity between the 82 nodes and subjected to analyses of node centrality and community structure. Whereas the consensus network was comprised of affective, integrative, and vestibular nodes, PCVD participants exhibited diminished integration and centrality among vestibular and affective nodes and increased centrality of visual, supplementary motor, and frontal and cingulate eye field nodes. Clinical outcomes, derived from dynamic posturography, were associated with approximately 62% of all connections but best predicted by amygdalar, prefrontal, and cingulate connectivity. No group-wise differences in diffusion metrics or tractography were noted. These findings indicate that cognitive, affective, and proprioceptive substrates contribute to vestibular processing and performance and highlight the need to consider these domains during clinical diagnosis and treatment planning.

INTRODUCTION

The term “vestibular” is applied to a variety of overlapping functional domains, including self‐orientation, sensory organization, sensory weighting, and static and dynamic balance. Previous studies suggest that the vestibular network includes the vestibular end organs, brainstem nuclei, and cortical areas including the multisensory orientation areas (MSO), located at the junction of the temporal and parietal cortices, and parieto‐insular vestibular area (PIVC), which is comprised of posterior insular complex (PIC) and parietal opercular area OP2 (Alsalman et al., 2016; Dieterich & Brandt, 2015; Indovina et al., 2020; Kirsch et al., 2018; Zu Eulenburg, Caspers, Roski, & Eickhoff, 2012). However, higher vestibular function, and some forms of vestibular dysfunction, may also involve feedback and feed‐forward cognitive processes as well as integration of multiple sensory modalities (Brandt, Strupp, & Dieterich, 2014). This implies that vestibular processing necessitates involvement of distributed brain regions subserving a variety of functions, which is supported by recent studies using probabilistic tractography. Raiser et al. (2020) and Indovina et al. (2020) report projections between known vestibular regions and frontal, prefrontal, parietal, cingulate, occipital, and posterior temporal cortices, as well as the medial temporal cortex and hippocampus.

Visual provocation of vestibular symptoms also elicits differential BOLD responses in these “ancillary” areas, which are, in turn, predictive of scores on clinical symptom scales. Riccelli et al. (2017), for example, demonstrated increased activity in the visual cortex with provocation in patients with persistent postural‐perceptual dizziness (PPPD), and noted that this recruitment was positively associated with scores on the Dizziness Handicap Inventory (Jacobson & Newman, 1990). An association between vestibular provocation and insular recruitment in controls, but not in PPPD participants, was reported in the same study (Riccelli et al., 2017). Similarly, it was recently reported that control participants and patients with subacute post‐concussive vestibular dysfunction (PCVD) exhibit differential recruitment of the visual, parieto‐insular, parietal, frontal, and cingulate cortices and hippocampus in response to visual‐vestibular provocation, and noted that activity in the frontal eye field, posterior hippocampus, and middle temporal visual area (hMST/V5) was associated with Vestibular Ocular Motor Screening (VOMS: Mucha et al., 2014) scores across all control and PCVD participants (Allen et al., 2021). Subjects with higher resting‐state connectivity between the posterior hippocampus and ipsilateral V5, and among the left and right‐hemisphere MSO, PIVC, and right anterior insula, have also been shown to exhibit increased symptoms during VOMS testing (Trofimova et al., 2021). These results not only implicate visuospatial and integrative areas in vestibular function and dysfunction, but also suggest that, in some cases, dysfunction may arise from an overreliance on visual cues. PCVD patients in particular appear to demonstrate greater recruitment of visual and visuospatial areas when provoked with visual stimuli containing high optic flow (Allen et al., 2021). They also exhibit richer, visuospatially directed networks, with “hubs” in the posterior hippocampus and anterior insula, and a greater number of intervening connections (longer path lengths) between the early visual areas V1/V2 and the posterior hippocampus (Trofimova et al., 2021).

Vestibular processing and dysfunction may also be intertwined with cognitive and affective states (Brandt et al., 2014). In addition to the induction or exacerbation of anxiety, depression, and other affective conditions co‐occurring with vestibular dysfunction, there is evidence that affective states themselves may influence the central vestibular systems. Staab (2014), for example, have proposed that state (momentary) and trait (continual) anxiety contribute to oculomotor gain and gaze control, and Passamonti, Riccelli, Lacquaniti, Staab, and Indovina (2018) have demonstrated that neuroticism scores in PPPD patients were positively associated with inferior frontal‐to‐extrastriate connectivity during visual‐vestibular provocation. Moreover, the insular, orbitofrontal, and cingulate cortices are also known to modulate affect (Boccia et al., 2016), and experimental manipulations of vestibular sensation have yielded consistent affective responses (Ponzo, Kirsch, Fotopoulou, & Jenkinson, 2018; Preuss, Hasler, & Mast, 2014). Vestibular functions may thus mediate body representation and ownership (Lopez, 2016; Mast, Preuss, Hartmann, & Grabherr, 2014) in addition to aiding spatial cognition and interoception.

It is therefore likely that more cortical regions contribute to vestibular processing and dysfunction than previously investigated. This may be similar to chronic pain, in which brain regions subserving the affective sequelae of pain are thought to mediate somatotopy in primary somatosensory cortex, as well as pain resilience, via first‐ and second‐order connections (Brooks, Zambreanu, Godinez, Craig, & Tracey, 2005; Flor, Nikolajsen, & Staehelin Jensen, 2006; Makin & Flor, 2020; Wei, Li, & Zhuo, 1999). These affective areas, as well as regional modulators of pain such as thalamus, cerebellum, and parietal cortex, have been grouped with somatosensory cortex as members of a “pain neuromatrix” (Hahamy & Makin, 2019; Jones & Pons, 1998; Kaas, Florence, & Jain, 1999; Zeharia, Hertz, Flash, & Amedi, 2015; Zeharia, Hofstetter, Flash, & Amedi, 2019). We posit that ancillary regions which form close connections with vestibular cortical areas may likewise be subsumed within a “vestibular neuromatrix” which may account for commonly observed clinical phenotypes characterized by persistent dizziness.

There were two objectives for the current study. First, we attempted to more completely define the “vestibular neuromatrix.” Prior studies of subjects with PCVD were based upon small sample sizes and relied on regions of interest (ROIs) either defined a priori from previously published work (Blum, Habeck, Steffener, Razlighi, & Stern, 2014; Della‐Justina et al., 2015; Eickhoff, Weiss, Amunts, Fink, & Zilles, 2006; Jamadar, Fielding, & Egan, 2013; Kolster, Peeters, & Orban, 2010; Poppenk, McIntosh, Craik, & Moscovitch, 2010; Vernet, Quentin, Chanes, Mitsumasu, & Valero‐Cabre, 2014; Zu Eulenburg et al., 2012) or derived a posteriori from differences in mean regional BOLD signal amplitudes between control and PCVD patients. Thus, here we expand the number and breadth of ROIs by incorporating (1) oculomotor areas, particularly those which subserve smooth pursuits and saccadic functions; (2) regions which exhibit differences in mean functional connectivity, rather than BOLD signal amplitude, between control and subacute PCVD subjects; and (3) regions which exhibit differences in mean functional connectivity as a function of personal emotional attributes or emotional regulation, including anxiety, depression, and post‐traumatic stress. Our second objective was to test the hypothesis that PCVD patients, who exhibit clinically significant scores on standard clinical instruments, including the Post‐Concussion Symptom Scale (PCSS: Lovell et al., 2006), Standardized Assessment of Concussion (SAC: McCrea, 2001; Borich et al., 2013), Balance Error Scoring System (BESS: Riemann, Guskiewicz, & Shields, 1999; Finnoff, Peterson, Hollman, & Smith, 2009; Bell, Guskiewicz, Clark, & Padua, 2011), and VOMS, are characterized by an overreliance on visual cues, and that PCVD therefore arises from pathological overweighting of visual information, which would be evinced by increased prominence of visual and multisensory integration areas within this neuromatrix.

METHODS

Study population and testing

This study was approved by our local Institutional Review Board, and all subjects provided informed consent. Then, 27 healthy control subjects (15 females, mean ± SD age 28.50 ± 6.44 years) and 18 patients with subacute PCVD (12 female, age 23.43 ± 4.90 years) were selected for participation in this study (see Table 1). Inclusion criteria for PCVD subjects included (a) a diagnosis of concussion, as defined by the World Health Organization Collaborating Centre for Neurotrauma Task Force (Carroll et al., 2004; Kristman et al., 2014), occurring 4–12 weeks prior to enrollment; and (b) clinical evidence for vestibular impairment, including subjective report of dizziness and/or imbalance or abnormal VOMS scores. Prospective participants with abnormal balance but no balance or VOMS abnormalities were not recruited. Exclusion criteria for all subjects were age <16 years or history of moderate or severe head injury, intracranial hemorrhage, seizure disorder, prior neurological surgery, peripheral neuropathy, musculoskeletal injuries affecting gait and balance, or chronic drug or alcohol use. Subjects with abnormal head impulse testing or videonystagmography consistent with peripheral vestibular hypofunction or benign paroxysmal positional vertigo were also excluded.

TABLE 1

Demographic and clinical data, presented as mean (SD) and median (intraquartile range). Normality was tested via Shapiro–Wilks test. If data were determined to follow a normal distribution, results are presented as “mean (SD)” and group differences assessed by independent‐samples t test, using unpooled sample variances if indicated by Levene's test. Otherwise, results are reported as “median (intraquartile range)” and group differences assessed by Mann–Whitney U test. Clinical data were missing from four control participants

	Control	PCVD	p‐Value
N (female)	27 (15)	18 (12)	n/a
Age	25.13 (2.36)	22.67 (3.35)	.003
Mean framewise displacement (mm)	0.12 (0.04)	0.10 (0.03)	.033
Maximum framewise displacement (mm)	0.51 (0.23)	0.68 (0.62)	.154
Lifetime concussions	0.00 (0.00)	2.00 (1.00)	n/a
Days since most recent concussion	n/a a	35.50 (29.00)	n/a
Clinical concussion measures
PCSS total	1.00 (4.00)	34.00 (27.00)	<.001
SAC	28.00 (3.00)	26.00 (3.00)	.003
BESS	10.00 (5.00)	12.00 (8.00)	.027
VOMS
Baseline	0.00 (0.00)	5.00 (3.00)	<.001
Smooth pursuits	0.00 (0.00)	5.00 (5.00)	<.001
Horizontal saccades	0.00 (0.00)	6.00 (6.00)	<.001
Vertical saccades	0.00 (0.00)	7.00 (7.00)	<.001
Near‐point convergence	0.00 (0.00)	6.00 (10.00)	<.001
Horizontal VOR	0.00 (1.00)	8.00 (7.00)	<.001
Vertical VOR	0.00 (0.00)	9.00 (10.00)	<.001
VMS	0.00 (1.00)	9.00 (9.00)	<.001
Affective measurements
Beck Depression Inventory	2.77 (4.34)	Not assessed b
Beck Anxiety Inventory	2.38 (3.86)	11.38 (8.68)	.002
PCL‐5	0.46 (1.127)	15.08 (7.46)	<.001
Objective vestibular/balance tests
SOT
Composite score	78.79 (5.30)	77.13 (6.44)	.141
Condition 2	92.25 (2.65)	89.25 (5.21)	.051
Condition 3	92.79 (2.59)	89.96 (3.96)	.043
Visual score	79.83 (10.08)	80.81 (11.05)	.468
Vestibular score	72.75 (13.27)	73.56 (9.63)	.504
Horizontal smooth pursuits
Mean velocity gain, 0.1 Hz	0.96 (0.91)	0.94 (0.89)	.202
Mean velocity gain, 0.5 Hz	1.00 (0.98)	1.00 (0.95)	.288
Mean phase, 0.1 Hz	2.19 (1.12)	1.46 (1.07)	.095
Mean phase, 0.5 Hz	0.96 (0.91)	0.94 (0.89)	.486
Velocity gain asymmetry, 0.1 Hz	2.11 (1.10)	2.73 (1.93)	.216
Velocity gain asymmetry, 0.5 Hz	1.33 (0.53)	1.85 (0.51)	.280
Horizontal random saccades
Peak velocity AUC	10,611.60 (1,239.50)	10,155.26 (1,037.60)	.669
Mean latency	0.19 (0.02)	0.19 (0.02)	1.000
Mean final accuracy	97.27 (6.93)	95.46 (5.38)	1.000

Note: Italics indicate data non‐normality as determined by Shapiro–Wilks test. Significant p values (<. 001) are represented in bold.

Abbreviations: BESS, Balance Error Scoring System; IQR, intraquartile range; PCL‐5, Post‐Traumatic Stress Disorder Checklist for the DSM‐5; PCSS, Post‐Concussive Symptom Scale; PCVD, post‐concussive vestibular dysfunction; SAC, Standardized Assessment of Concussion; SD, standard deviation; SOT, sensory organization test; VMS, visual motion sensitivity; VOMS, Vestibular Ocular Motor Screening; VOR, vestibulo‐ocular reflex.

One control participant reported a concussion which occurred 19 years prior to the study.

Beck Depression data were obtained from an insufficient number of PCVD participants for evaluation of group differences.

Clinical assessments, including PCSS and VOMS, were used to assess vestibular symptoms. In addition to the VOMS, balance and vestibular assessment included the BESS and the objective sensory organization test (SOT) obtained via dynamic posturography in an immersive virtual projection system (Bertec FIT CDP/IVR, Bertec Immersive Labs, Columbus, OH) equipped with dual‐balance force plates. SOT has established validity for vestibular impairment after stroke, traumatic brain injury, and vestibular disorders (Grove, Whitney, Hetzel, Heiderscheit, & Mark Pyle, 2021; Nashner, 1993). Briefly, the SOT is comprised of six sensory conditions which test a participant's ability to maintain equilibrium, operationalized by postural sway and center of pressure, over three 20‐s trials per condition (Monsell, Furman, Herdman, Konrad, & Shepard, 1997; Nashner, 1993). Conditions include (a) eyes open while on a fixed support surface and facing a static immersive scene; (b) eyes closed on a fixed surface; (c) eyes open with sway‐referenced visual surround (i.e., the visual scene changes synchronously with variations in center of gravity); (d) eyes open on a sway‐referenced support surface (i.e., the support surface changes synchronously with variations in center of gravity); (e) eyes closed on a sway‐referenced support surface; and (f) eyes open on both a sway‐referenced support surface and visual surround. For each of the six conditions, an equilibrium score is generated from the center of gravity sway and an overall SOT composite score is reported representing performance on all six conditions, with higher scores representing better balance performance. In addition, standardized sensory subscores are calculated representing utilization of specific visual, vestibular, or somatosensory sensory strategies. Group‐wise differences in SOT scores were assessed by Mann–Whitney U test.

In addition, we administered the SAC, a bedside cognitive screen, and three self‐report affective assessments, including the Beck Depression Inventory, second edition (BDI‐II: Beck, Steer, & Brown, 1996; Dozois, Dobson, & Ahnberg, 1998; Richter, Werner, Heerlein, Kraus, & Sauer, 1998), Beck Anxiety Inventory (BAI: Steer & Beck, 1997; Ulusoy, Sahin, & Erkmen, 1998), and Post‐traumatic Stress Disorder Checklist for the DSM‐5 (PCL‐5: Blevins, Weathers, Davis, Witte, & Domino, 2015; Zuromski et al., 2019). Group‐wise differences on demographic and clinical measures were evaluated by t test or Mann–Whitney U, if the Shapiro‐Wilks test indicated a violation of normality (see Table 1).

Magnetic resonance imaging acquisition

Magnetic resonance imaging (MRI) was conducted on a 3 T Siemens Prismafit MRI with 32‐channel head coil. −weighted anatomical data were obtained as 3D‐multiecho MPRAGE with inversion recovery (van der Kouwe, Benner, Salat, & Fischl, 2008); sagittal slices; 1.0 mm isotropic resolution; 256 mm field of view (FOV); TR/TE 2,300/2.96 ms, 7° flip angle; 1,100 ms inversion). ‐weighted resting‐state functional MRI (rsfMRI) data were acquired as simultaneous multislice echo planar images (Feinberg et al., 2010; Moeller et al., 2010; Xu et al., 2013) over a continuous 7 min 30 s acquisition (600 volumes, multiband factor 5, 2.5 mm isotropic resolution with 50 axial slices, anterior‐to‐posterior phase encoding, 220 mm FOV, TR/TE 750/32 ms, 52° flip angle, matrix size = 88 × 88). For diffusion imaging and tractography, we obtained multishell diffusion weighted imaging with optimal angular coverage using 128 diffusion directions distributed over four shells (4 volumes of ; 17 of ; 39 of ; and 68 volumes of ) with 2 mm3 isotropic voxel resolution; multiband factor 3; 79 ms TE; 2,750 ms TR; 78° flip angle; AP phase encode; 232 × 256 FOV. Additionally, 12 images were acquired interspersed between the diffusion volumes. We also acquired two volumes of in the opposite phase encoding direction to correct for susceptibility‐based distortion using the FSL TOPUP utility.

Diffusion data preprocessing and analysis

Our preprocessing for diffusion modeling and tractography (described below) consisted of standard workflows in the FMRIB Software Library (FSL version 6.0: Wellcome Centre for Integrative Neuroimaging, Oxford, UK) which have been described in detail elsewhere (Behrens et al., 2003; Behrens, Berg, Jbabdi, Rushworth, & Woolrich, 2007; Jbabdi, Sotiropoulos, Savio, Grana, & Behrens, 2012). In brief, these steps consisted of fieldmap‐facilitated distortion correction via TOPUP, diffusion tensor model fitting in DTIFIT, and Bayesian estimation of diffusion parameters, including modeling of crossing fibers, in BedpostX. Additionally, we projected subject data onto a mean tractography skeleton using TBSS (Smith et al., 2006), and assessed voxelwise differences in fractional anisotropy and mean diffusivity between control and PCVD subjects via two‐sample unpaired t tests in FSL randomize (Winkler, Ridgway, Webster, Smith, & Nichols, 2014). Cluster‐level inferences on these results were established using threshold‐free cluster enhancement (TFCE; Smith & Nichols, 2009).

rsfMRI preprocessing

Preprocessing and bivariate correlation (seed‐to‐voxel and ROI‐to‐ROI connectivity) analyses were performed in the CONN Toolbox v19c (nitrc.org/projects/conn) (Whitfield‐Gabrieli & Nieto‐Castanon, 2012), which wraps SPM8 (www.fil.ion.ucl.ac.uk/spm) and aCompCor (Behzadi, Restom, Liau, & Liu, 2007) noise source removal functions. In addition to standard slice timing, field map, and motion correction preprocessing, anatomical and rsfMRI datasets were registered to the Montreal Neurological Institute (MNI152 stereotactic space) template via SPM's nonlinear, diffeomorphic method (DARTEL), which afforded both forward and backward transformations between subject and MNI template space (Ashburner, 2007; Ashburner & Friston, 2011). A priori MNI‐space tissue probability maps were then applied to the data for initial classification of voxel tissue types as gray matter, white matter, or CSF, and the tissue estimates refined by direct segmentation of the anatomical image by Gaussian mixture model. Finally, as part of the aCompCor workflow, the motion parameters, their first‐order derivatives, ART‐defined TR outliers and mean CSF‐ and white matter‐derived signals were regressed from each subject's rsfMRI dataset via ordinary least squares (OLS) regression (see Whitfield‐Gabrieli & Nieto‐Castanon, 2012 for additional details). We then smoothed the residual data with an 8 mm FWHM Gaussian kernel and, in consideration of our 750 ms sampling rate and to avoid filtering potentially meaningful higher‐frequency content (Boubela et al., 2013), applied band‐pass filtering at 10–250 MHz (Birn, Diamond, Smith, & Bandettini, 2006; Fox, Snyder, McAvoy, Barch, & Raichle, 2005; Fox, Zhang, Snyder, & Raichle, 2009). All structural and “denoised” functional data, gray matter, white matter, and CSF masks, as well as ROIs (see below), were manually inspected to confirm registration validity. We also present the outputs of our quality control workflow as Supplementary Figure S1. Group‐wise mean and maximal motion (framewise displacements) are also reported in Table 1. Maximal framewise displacement, across all subjects, was 0.76 mm, and a maximum nine TRs excluded from any one subject by the aCompCor workflow.

Network node selection

We selected 82 ROIs for this study based on previous publications (Berman et al., 1999; Blum et al., 2014; Della‐Justina et al., 2015; DeSouza, Menon, & Everling, 2003; Dumoulin et al., 2000; Eickhoff et al., 2006; Jamadar et al., 2013; Kellar, Newman, Pestilli, Cheng, & Port, 2018; Kolster et al., 2010; Poppenk et al., 2010; Vernet et al., 2014; Zu Eulenburg et al., 2012) and empirical findings from the current dataset (see Figure 1 and Supplementary Table ST2). Twenty‐five of these ROIs corresponded to putative cortical substrates of vestibular, spatial localization and orientation, and saccadic functions, and were defined in the MNI stereotaxic space as 8‐mm‐diameter spheres centered at coordinates based on coordinates from existing literature. Although most of these sources provided ROI coordinates in MNI space, a minority did so in Talairach–Tourneaux stereotactic space, and were converted to MNI using the MNI2TAL service (Yale University: https://bioimagesuiteweb.github.io/webapp). Three ROIs, which we localized to the right anterior insula and visual areas V1/V2 and V5, were included based on differences in BOLD signal magnitude between control and PCVD participants, as indicated by t test, as they did not overlap with the 25 a priori ROIs.

FIGURE 1

Source regions of interest (regions of interest [ROIs] or nodes), derived a priori from previous literature, red, pink; t tests, orange, or multivoxel pattern analysis (MVPA) on group, green; or MVPA against Beck Anxiety Inventory (BAI, cyan), Beck Depression Inventory (BDI, blue), or post‐traumatic stress metrics (PCL5, violet)

The remaining ROIs were derived from the current dataset via multivoxel pattern analyses (MVPAs) to assess the whole‐brain, multivariate pattern of pairwise voxel connections and to localize voxel clusters which differed in their average whole‐brain connectivity (rather than BOLD signal magnitude, as described above) as a function of group (control vs. PCVD; 22 ROIs). Although several of the ROIs derived from MVPA on group, including medial prefrontal, anterior cingulate, insular, orbitofrontal cortex, and the hippocampus, are associated with depression (Gudayol‐Ferré, Peró‐Cebollero, González‐Garrido, & Guàrdia‐Olmos, 2015) or anxiety (Mochcovitch, da Rocha Freire, Garcia, & Nardi, 2014), additional affective ROIs were derived for our current subject group via separate MVPAs against participants' scores on depression (BDI, 4 ROIs), anxiety (BAI, 21 ROIs), or other aspects of emotional regulation (PCL‐5 checklist for the DSM‐5 [PCL5], 16 ROIs). In each case, MVPA analyses were conducted in the CONN Toolbox, which implements MVPA (Beaty, Silvia, & Benedek, 2017; Hankin, Badanes, Smolen, & Young, 2015; Whitfield‐Gabrieli et al., 2016; Whitfield‐Gabrieli & Nieto‐Castanon, 2012) by first computing each subject's cross‐covariance matrix C across all gray matter voxels, then performing a singular value decomposition on C to produce a voxel‐to‐voxel correlation matrix R, such that the norm of one row of R represents the strength of connectivity, for one participant, between a given voxel and all other voxels in the brain. The dimensionality of this multivoxel pattern across all participants is then reduced via principal components analysis, and the top n principal components selected for general linear modeling, where n is generally taken to be a 5:1 participant/component ratio, that is, six components for the current study (Hair, Black, Babin, Anderson, & Tatham, 2006). Finally, the selected principal components are subjected to a general linear model to obtain an ‐statistic for each effect of interest among group, BDI, BAI, and PCL5. For this exploratory omnibus analysis, and in anticipation of more stringent statistical thresholds in post hoc analyses, we opted to conduct these GLMs with liberal cluster growth and peak voxel thresholds of and (uncorrected), respectively. ROIs were defined in MNI space as 8 mm‐diameter spheres centered at the peak voxel of the suprathreshold clusters; those which overlapped with other a priori or a posteriori ROIs, and those which were localized to, or possessed a peak voxel in, non‐gray matter areas were excluded. These ROIs were back‐projected to subject space and manually inspected for registration validity. We then extracted the mean ROI timeseries from the denoised but unsmoothed rsfMRI data of each subject and subjected them to OLS regression to obtain a bivariate ROI–ROI correlation coefficient. Connectivity values were Fisher‐transformed to Z‐scores and stored in a functional connectivity matrix (FCM: dimensions ) for each subject.

Computation of “consensus” and group‐specific networks

Generation of functional connectivity and functional adjacency matrices

Graph‐theoretic analyses are often predicated upon sparse or binary‐valued adjacency matrices, in which graph edges—ROI–ROI connections, in this case—are assigned a Boolean truth‐value depending on whether or not their weights (connection strengths) exceed a predetermined threshold. We therefore obtained “functional adjacency matrices” (fAMs) by a workflow consisting of four steps. First, we conducted one‐sample t tests on subjects' FCMs using the MATLAB function ttest to obtain uncorrected p‐values. Next, we performed a series of simulations (see Supplementary Figure S3) to establish an acceptable connection‐level statistical threshold while correcting for multiple comparisons by Benjamini–Hochberg false discovery rate, or (Benjamini & Hochberg, 1995) or Storey's bootstrapping method, or (Storey, 2002). Based on the number of surviving connections as a function of corrected p‐value in these simulations (Supplementary Figure S3), the lenient thresholds used to derive some of the ROIs using MVPA, the exploratory nature of this study, its small sample size, and a need to minimize connection‐wise false‐positive rate, we selected a conservative threshold of for further analysis. We then conducted separate one‐sample t tests on control (HC) and PCVD (ST) subjects and applied the criterion to each group's t‐test results. Finally, we produced two‐dimensional, binary‐valued fAMs for each group (HC fAM and ST fAM) comprised of “1” entries for each suprathreshold connection and “0” entries for all other connections (Figure 2a). All further analyses were conditioned on the HC and ST fAMs. The FCMs and fAMs are presented in further detail as Supplementary Figure S4.

FIGURE 2

Generation of the structural (sAM) and functional (fAM) adjacency matrices (a,b); control (HC) and post‐concussive vestibular dysfunction (PCVD) (ST) “core” networks (c); and consensus network (d). See text for the details of each step

Tractography and generation of structural adjacency matrices

In order to more rigorously constrain the group‐wise fAMs with structural data by excluding functional “connections” which were insufficiently supported by physical connections, we leveraged fMRI‐guided, region‐local probabilistic tractography (reviewed in (Zhu et al., 2014)) using probtrackx2 (Behrens et al., 2007) to estimate anatomical connectivity among the 82 ROIs. Then, 5,000 streamlines were initiated from each seed voxel using each of the remaining 81 ROIs as termination masks, yielding a structural connectivity matrix for each subject. Each cell in these matrices represented the number of successful replicates of a given connection out of these 5,000 iterations. No group‐wise differences in fractional anisotropy or mean diffusivity were observed after TFCE, , and permutation tests (permutationTest.m MATLAB code: Laurens Krol, Berlin Institute of Technology) failed to detect any differences in structural connectivity between control and PCVD subjects, , after 5,000 iterations, indicating that HC and ST participants did not differ with respect to structural connectivity. Consequently, we computed the mean structural connectivity over all subjects and iteratively increased the minimal number of replicates (the “threshold”: see Supplementary Figure S3) until obtaining a mean nodal degree = 12, based on findings by van den Heuvel and Sporns (2011). This “threshold” rejected the lowest 3.10% of replicates (<155 replicates out of a maximum 5,000) and yielded a network density of approximately 20% (Rosen & Halgren, 2021). Connections which exceeded this replication threshold were assigned a value of in an overall structural adjacency matrix (sAM: Figure 2b); all others received a value of . The SCM and sAM are presented in further detail as Supplementary Figure S4. Alternative thresholding of the SCM at 70% sparsity, or retaining the top 30% of connections (Buchanan et al., 2020) yielded comparable results, although some structural connections were retained by the nodal degree, but not sparsity, method, suggesting that connections with the highest number of replicates were generally stable regardless of the threshold method used (see Supplementary Figures S3 and S14).

Fusion of structural and functional connectivity matrices

As a final preparatory step, we leveraged the all‐subjects structural and group‐wise fAMs to calculate a sparse overall “consensus” matrix as well as group‐specific sparse “core” matrices (the “HC core” and “ST core”). The consensus matrix included connections common to both control (HC) and PCVD (ST) subjects, whereas the HC core included connections from the consensus matrix as well as connections which were suprathreshold in the HC subjects but not ST, and likewise for the ST core. We computed the HC and ST cores as a Boolean union of the HC or ST fAMs and the sAM, that is, [HC fAM sAM] and [ST fAM sAM], respectively, where represents a logical AND operation, equivalent to elementwise multiplication by 0 or 1 (Figure 2). The consensus adjacency matrix was comprised of connections common to both HC and ST subjects which also survived masking by the sAM, that is, [HC fAM ST fAM sAM] (Figure 2d). Note that, because the SAM was derived from the mean SCM computed across all participants, the same SAM was applied to both HC and ST fAMs. Differences in core networks were therefore due to differences in the fAM component.

Statistical and graph‐theoretic assessments of subacute PCVD

We leveraged a combination of graph‐theoretic and mass‐univariate tests to (a) isolate group‐wise differences in connectivity strengths within the consensus network, (b) determine the degree of integration of visual and ancillary ROIs into the vestibular network, and (c) test our hypothesis that PCVD emerges from overweighting of visuospatial substrates. Our graph‐theoretic measures of interest for the present study included betweenness centrality (consensus and core networks) and nodal community structure (within the core networks) and were intended to address objectives (b) and (c). The betweenness centrality is a measure of the “hub‐ness” of a given node (ROI) in a graph; more rigorously, it is defined as the number of times a node appears on the shortest paths between two nodes and . As our analysis was performed on undirected, unweighted graphs derived from the consensus and core adjacency matrices, the “shortest path” was computed by breadth‐first algorithm. Noting that there may be multiple “shortest paths” between and , we take the average of the number of paths between and which contain node [], over the total number of shortest paths between and ():

These nodal centralities were normalized to the total number of nodes in the graph, :

Nodal community structure was determined for the HC and ST core networks via Louvain algorithm (community_louvain utility in the Brain Connectivity Toolbox [BCT]: brain-connectivity-toolbox.net (Rubinov & Sporns, 2010)), which subdivides a graph into topological neighborhoods (“cliques” or “communities”) by iterative maximization of the number of edges within a clique and minimization of the number of edges between cliques. The Louvain implementation in the BCT is deterministic with a default resolution parameter () of 1. Solutions with alternative values from 0 (more weighting toward smaller modular communities) to 10 (weighting toward larger modular communities), averaged over 5,000 iterations for each , are provided as Supplementary Figure S15. Full descriptions and derivations of this algorithm are provided in Reichardt and Bornholdt (2006) (Reichardt & Bornholdt, 2006) and Rubinov and Sporns (2011) (Rubinov & Sporns, 2011).

To address objective (a)—the assessment of group differences in connection strengths within the consensus network—we masked the HC and ST FCMs by the consensus adjacency matrix (i.e., [HC FCM HC fAM ST fAM sAM] and [ST FCM HC fAM ST fAM sAM]) and conducted independent‐samples t tests on the extant connections. These t tests were corrected for multiple comparisons by Storey's bootstrap method at .

Regression of network weights against clinical measures

Associations between SOT condition or composite scores and functional connectivity within the core networks were assessed via OLS regressions using the Python package statsmodels (v. 0.12.2, www.statsmodels.org). Separate regressions were conducted for each connection in the control or PCVD core networks, each of which included consensus (shared) connections. Inclusion of the covariates days elapsed since injury and age as covariates incurred rank deficiency and we observed that the regression coefficients for these terms were one to two orders of magnitude less than those for connectivity or group. Consequently, the final regression models included only group (encoded as a binary categorical variable, HC = 0, ST = 1) and the normalized bivariate correlation coefficient for the connection as predictors as well as an intercept term.

RESULTS

Demographics and clinical scores

Demographic and clinical scores are shown in Table 1. PCVD participants (; 12 female, median age 23 years) presented with a median of two lifetime concussions and 35.5 days elapsed between the most recent concussion and MRI scan. Conversely, of our control participants (; 15 females, median age 25 years), only one noted a previous concussion, which occurred 19 years previously. Group and sex differences with respect to age and number of lifetime concussions were evaluated by general linear model and indicated that, although control and PCVD groups differed with respect to age , , , participant age was not associated with the number of lifetime concussions , , , . Mann–Whitney rank sum tests failed to resolve group‐wise differences in SOT conditions 1, 4, 5, or 6, SOT vestibular subscore, or visual subscore, but indicated that PCVD participants did exhibit greater sway variability in conditions 2 (eyes‐closed, fixed surface; , ) and 3 (sway‐referenced visual surround; , ). Three control and two PCVD participants did not complete posturography assessments and were not included in the regressions against core network connectivity described below.

Consensus and “core” networks

The consensus network is depicted in Figure 3. In addition to left and right MSO, PIVC, and right supramarginal gyrus (SMG), the network consisted of cingulate (anterior, middle, and subgenual), prefrontal and orbitofrontal areas, the frontal and supplementary eye fields, motor, premotor, and supplementary motor cortices; calcarine, lingual, and dorsal visual areas (inferior parietal sulcus ROIs believed to subserve eye movements, the superior parietal lobule [SPL], visual area V5 in the medial temporal region and V6 in the parieto‐occipital sulcus region), fusiform cortex, angular gyrus, and cerebellar lobule VIIa. Anterior and subgenual cingulate and prefrontal Brodmann area (BA) 47 formed consensus connections with temporal and mesial temporal areas including the left amygdala and rostral and posterior hippocampi. Prefrontal BA 47, the dorsal aspect of anterior cingulate area 32 (R.ACC.32d), the dorsal aspect of middle cingulate BA 23 (L.MiCC.23d), ventromedial BA5, and the medial aspect of right superior parietal BA 7 (R.SPL.7 m) exhibited high centralities in the consensus network. The “core” network of control participants included functional connections between right anterior insula and the ipsilateral MSO and SMG, as well as between prefrontal BA 45, left MSO, and middle cingulate, which did not survive thresholding in PCVD subjects. Additionally, functional connectivity between the right caudate nucleus and ipsilateral anterior cingulate, as well as between the left caudate nucleus and ipsilateral posterior hippocampus and amygdala, did not survive thresholding in the PCVD group. No group‐wise differences in structural connectivity or diffusion metrics, however, were detected. We also noted that the centralities or “hub‐ness” of the middle and anterior cingulate, left MSO, and right‐hemisphere anterior insular, inferior parietal, and orbitofrontal cortices were markedly diminished in the PCVD group. Conversely, the centrality of primary and secondary visual and supplementary motor cortices, posterior (but not rostral) hippocampus, and the frontal eye fields and supplementary cingulate eye areas (SCEA) were increased in the ST group, as shown in Supplementary Figure S9. Nodal community structures, which are presented as Supplementary Figures S5 and S6, indicated that vestibular and cognitive–affective ROIs were closely associated, whereas disjunctions between HC and ST communities (Supplementary Figure S7) suggested group‐wise differences in the integration of cingulate, sensorimotor, frontal eye field, and posterior hippocampal areas, among others, into functional networks. (Note that community membership may differ between groups even for “consensus” connections; this simply indicates that the connection is not as tightly integrated into a topological neighborhood in one group as it is in another.) These structures were derived using a Louvain resolution parameter of , which gives equal preference to small and large communities. However, nodal community assignments were observed to vary with various values from 0 (weighting toward smaller communities) to 10 (weighting toward larger communities), such that assignments for a given node tended to become less stable, and overall modularity (Q‐metric) tended to decrease, with larger values. The greatest variation in community assignments was observed for sensorimotor, early visual, prefrontal, and parietal and parieto‐occipital nodes. These results are presented as Supplementary Figure S15.

FIGURE 3

Network edges common to all subjects (“consensus network,” Panel a) versus the networks of control (“HC,” Panel b) and post‐concussive vestibular dysfunction (PCVD) participants (“ST,” Panel c). Node size in all panels corresponds to betweenness centrality. Yellow edges in (b) and (c) denote consensus edges; violet edges in (b) were suprathreshold in HC, but not ST, and in (c), suprathreshold in ST, but not HC. Anterior and middle cingulate (ACC, MiCC) and the precuneus (PCu) demonstrate high betweenness centrality across subjects; however, note the increased centrality of cerebellar, early visual, and dorsal visual stream areas in ST (c). Node labels are provided as Table ST2 in the Supplementary Materials

Permutation tests on the core network connections indicated several connections which were stronger in the control group than in PCVD (), including left middle frontal gyrus—right frontal eye field (L.MiFG:R.FEF), right frontal eye field—right prefrontal cortex (ventral Brodmann BA 9/rostral 46/47), and right inferior frontal gyrus—right inferior temporal area TE1 (posterior). The PCVD group, however, exhibited stronger functional connectivity among orbitofrontal, prefrontal, and cingulate ROIs in the right hemisphere (rostral BA 11 and posterior BA 10, anterior BA9/ventral BA46, and dorsal BA 32 in the anterior cingulate); among cingulate, temporal, and medial temporal ROIs in the left hemisphere (posterior BA 33 in the anterior cingulate, posterior BA32 in the middle cingulate, putamen, amygdala, and left anterior BA38 in the superior temporal gyrus); among occipital and parietal structures in the right hemisphere (dorsomedial parieto‐occipital sulcus/area V6, angular gyrus, area V5, and medioventral BA 37 in the occipital fusiform gyrus); among the supplementary eye fields, supplementary motor cortex, and cingulate structures along the midline; and between right‐hemisphere MSO and the ipsilateral SMG. Of these, group‐wise differences in the right inferior frontal—right TE1p (HC > ST), right frontal eye field—right prefrontal ventral BA9/rostral 46/47 (HC > ST), left supplementary eye field—middle cingulate (ST > HC), and right MSO—right SMG connections attained the greatest statistical significance (Figure 4). We note some degree of discrepancy between these permutation test results and those of the one‐sample t tests conducted during the generation of the consensus and core adjacency matrices. As shown in Figure 5, positive skew in the distribution of connectivities within the PCVD group effected a rightward shift in the group mean. Consequently, many of the connection‐level permutation tests failed to refute the null hypothesis, whereas the one‐sample adjacency matrix tests did not.

FIGURE 4

Results of permutation tests on connection strength between control subjects (HC) and subjects with post‐concussive vestibular dysfunction (PCVD), corrected for multiple comparisons at . The ‐values for the analysis are presented as Panel (a). Connections which were stronger in HC than ST are presented in Panel (b); those which were stronger in ST are presented in Panel (c). Color bars refer to edge colors; node colors are as in Figure 3. Node labels are provided as Table ST2 in the Supplementary Materials

FIGURE 5

Demonstration of a “core” connection (left posterior hippocampus―left amygdala) in the control group which was not present in the “consensus” network or “core” network of the post‐concussive vestibular dysfunction (PCVD) group, and for which the group‐wise t test on connectivity was not significant. In this case, the distribution of PCVD participants' connectivities exhibits positive (leftward) skew, with a longer right tail. This moves the PCVD distribution's mean toward the right, and the group‐wise t test on connectivity was not significant. Conversely, subjecting each group's connectivities to a one‐sample t test against a mean of zero, as was performed during the creation of the adjacency matrices, indicates that the connection was “present”—the normalized bivariate correlation of the L.HCp and L.AMG signals differed significantly from zero—in the control group, but not in PCVD. HC, controls; ST, PCVD subjects; L HCp, left posterior hippocampus; L AMG, left amygdala

Regression results versus SOT scores

Detailed results of our OLS regression analyses are provided as Supplementary Figures S10–S12 and presented in summary form as Supplementary Table ST2. Here, we highlight the results for connections which yielded results with high statistical significance and, noting that we failed to discern group‐wise differences in SOT composite scores or in scores for conditions 1 or 4–6, will focus on regressions of consensus connections against SOT conditions 2 and 3.

Here, 65 of the 105 consensus connections (61.90%) were predictive of SOT condition 2 (eyes‐closed, fixed surface) score, condition 3 (sway‐referenced vision) score, or both, when leveraged in conjunction with the categorical variable group. Then, 33 of these 65 connections were significant predictors of both scores with an average values of .172 and .177 for conditions 2 and 3, respectively; 63 of the 65 were predictors of condition 3 scores (average ), but not condition 2 scores; and two connections (R.PFC.9a46v:R.PFC.9v46r47 and R.FEF:R.OFC.10p) were significant predictors of condition 2, but not condition 3, scores (average .147). The best model for condition 2 () consisted of the L.MiFG:R.IFG connection, which was comprised of a priori vestibular ROIs previous characterized as subserving visual/vestibulomotor integration (Della‐Justina et al., 2015), and was associated with a 0.89‐point decrease in SOT condition 2 score with every 0.10‐unit increase in functional connectivity between the left middle frontal gyrus and right inferior frontal gyrus (connection , group , , , ). This connection was also predictive of condition 3 score, but less significantly so (connection , group , , , ). The two best models for condition 3 score consisted of the connections L.FEF:L.PFC.45r (connection , group , , , ) and L.IFJa:L.MiFG (connection , group , , , ). These nodes were comprised of vestibular ROIs (L.FEF and L.MiFG: see (Berman et al., 1999; Della‐Justina et al., 2015; DeSouza et al., 2003; Jamadar et al., 2013; Kellar et al., 2018; Vernet et al., 2014)) as well as “affective” ROIs (PFC.45r, derived from an MVPA on PCL5 score, and IFJa, the anterior aspect of the inferior frontal junction, derived from an MVPA on BAI score). In both cases, an increase in connection strength was associated with a decrease in SOT condition 3 score. Both connections were also successful predictors of condition 2 score, but less significantly so (L.FEF:L.PFC.45r: connection , group , , , ; L.IFJa:L.MiFG: connection , group , , , ). This inverse relationship between functional connectivity strength and SOT condition score was generally persistent across all consensus connections, and strongest in the association between SOT condition 2 score and the L.AMG:L.SgCC.32 connection (left amygdala―left subgenual cortex, connection , model ), SOT condition 2 score and the R.PCC.31pd:L.MiCC.23d connection (right posterior cingulate―left middle cingulate, connection , model ), and between SOT condition 3 score and the L.FEF:L.PFC.45r connection (described above). Exceptions to this observation included connections among the anterior cingulate, middle cingulate, SCEA, prefrontal cortex, right caudate and putamen, left supplementary eye fields, left hippocampus, right angular gyrus, and right primary and supplementary motor areas. These connections were associated with 0.21‐ to 0.55‐point decreases in SOT score with every 0.10‐unit increase in connectivity strength.

DISCUSSION

There is considerable interest and ongoing effort to characterize the cortical network subserving vestibular function and recent clinical and neuroimaging studies evince complex interactions between central vestibular substrates and cognitive, oculomotor, and emotional/affective domains (Indovina et al., 2020; Raiser et al., 2020). Brain regions involved in visuospatial, affective, and integrative functions also subserve vestibular processing, modulate oculomotor gain and visual attention, and may contribute to dysfunction in a subset of syndromes including PPPD and PCVD (Allen et al., 2021; Passamonti et al., 2018; Riccelli et al., 2017; Trofimova et al., 2021). This study is, to our knowledge, the first to situate these ancillary connections in the larger context of vestibular function and dysfunction. The findings presented here represent additional structural and functional evidence that such “ancillary” areas are functionally integrated with vestibular regions and that many of these connections are predictive of vestibular performance (as measured by the SOT).

We have delineated a “consensus” network which not only includes vestibular areas such as posterior insula (PIVC) and multimodal integration areas (MSO, angular gyrus) but also regions subserving visuospatial processing, attention, and gain control (primary visual, dorsal visual, FEF, SPL/precuneus, supplementary cingulate eye fields, and fusiform cortex), as well as motor (primary, supplementary, and premotor cortex and cerebellum) and affective processes (anterior, subgenual, and middle cingulate, orbitofrontal cortex). Of these, orbitofrontal BA 47, ventromedial postcentral area 5, dorsal anterior and dorsal middle cingulate, and superior parietal area 7 exhibited high betweenness centrality, indicating that these regions serve as “hubs” of this consensus network. By definition, this network is common to both control participants and those with PCVD, and we therefore propose that this network may represent the necessary substrate for higher vestibular function. Despite methodological differences, several connections reported in Raiser et al. (Raiser et al., 2020), such as inferior parietal sulcus (“L.IPS.SP” and “R.IPS.SP” in the present study, vs. “IPS3” in Raiser), PIVC, SPL (“PCu” vs. “Area 7”), rostral/inferior precuneus (“L.SPL.7r5” and “R.SPL.7 m” vs. “VPS”), medial superior temporal areas (“V5.MT” and “MSO” vs. “hMST” and “PFcm.VPS”), posterior cingulate/visual cingulate cortex (“PCC” vs. “CSv”), and premotor and supplementary motor cortex, were centrally incorporated in the proposed consensus network. Several of these areas are known to subserve vestibular or saccadic functions (Berman et al., 1999; Kellar et al., 2018), and the nodal community structure reported by Raiser and colleagues comports well with our own findings. Due to our inclusion of “ancillary” ROIs, our present network is a superset of those previously reported. These ancillary ROIs appear to be closely integrated into nodal communities of vestibular and saccadic substrates (see Supplementary Figures S5–S7). We must note, however, that these community structures are contingent upon the specification of the Louvain resolution parameter , and that affiliations between nodes varied with different values. Somatomotor, early visual, parietal, and parieto‐occipital nodes were particularly susceptible to variation in community assignment as a function of (Supplementary Figure S15). Consequently, further replicative work is required to determine whether the integration of “vestibular” and “ancillary” nodes is generalizable across studies.

If the consensus network represents the fundamental substrate of higher vestibular function, the control and PCVD “core” networks may represent differential modulation of vestibular processing and output between “baseline” and “dysfunctional” states upon provocation. A comparison of network centralities between these two states suggests a shift from a prefrontal/cingulate/precuneus processing stream in controls to a superior parietal/visual/supplementary motor and cerebellar processing stream in PCVD. Control participants exhibited tight functional integration (in the “HC core” network) between the right‐hemisphere anterior insula, SMG, and MSO and between left‐hemisphere lateral prefrontal area 45, MSO, and middle cingulate cortex, whereas the connectivity of the insula―supramarginal―MSO motif was markedly diminished, and lateral prefrontal―middle cingulate―MSO motif largely absent, in PCVD participants. PCVD participants also exhibited (in the “ST core” network) decreased posterior hippocampus―amygdala and anterior cingulate―caudate connectivity and decreased centrality with respect to the dorsal anterior and dorsal middle cingulate, left MSO, and orbitofrontal “hubs,” which were partially supplanted by markedly increased centrality in primary and secondary visual and supplementary motor cortices, frontal eye fields, SCEA, and posterior hippocampus. Notably, differences in performance on SOT conditions 2 (eyes closed, fixed surface) and 3 (sway‐referenced vision), taken across all participants, were associated with differences in the connectivity of the amygdala, lateral prefrontal area 45, and subgenual and middle cingulate cortex, among others (Supplementary Table ST2).

These findings indicate that vestibular dysfunction may be associated with differential modulation of attention (Eckert et al., 2009; Touroutoglou, Hollenbeck, Dickerson, & Feldman Barrett, 2012), arousal (Rudebeck et al., 2014), proprioceptive (Ben‐Shabat, Matyas, Pell, Brodtmann, & Carey, 2015; Kheradmand, Lasker, & Zee, 2015) and somatomotor functions as well as movement planning (Oane et al., 2020; Potok, Maskiewicz, Kroliczak, & Marangon, 2019). Several of these substrates, including subgenual cortex (Hajek, Kozeny, Kopecek, Alda, & Hoschl, 2008; Hamani et al., 2011), right anterior insula (Hatton et al., 2012), SMG (Silani, Lamm, Ruff, & Singer, 2013), and anterior cingulate (Baiano et al., 2007; Bush, Luu, & Posner, 2000) are associated with affective and social cognitive processing in addition to attention, arousal, and somatomotor functions. These dual roles may more abstractly indicate functions in egocentric orientation and motor planning through body representation and ownership (Lopez, 2016; Mast et al., 2014) that may be disrupted in vestibular syndromes and drive clinical symptoms. More generally, however, our findings emphasize not only the complex nature of such syndromes but also the clinical relevance of their cognitive, affective, and proprioceptive sequelae. As behavioral performance within these domains may determine clinical outcomes, they must be systematically evaluated during clinical diagnosis and considered during treatment planning.

There are, however, several limitations in our study which should be addressed with further research. First, despite the expanded corpus of ROIs in the present study, we expect that we have failed to include other potentially important nodes in our proposed networks: for example, brainstem or thalamic areas were not included in the current work despite their known associations with vestibular processing, nor was a left‐hemisphere PIVC homolog ROI included. Our conclusions should therefore be viewed as an argument that the networks derived in this article participate in a “vestibular neuromatrix,” but likely do not represent a comprehensive description. Second, whereas our proposed neuromatrix is contingent upon functional connectivity metrics derived from normalized cross‐correlations between the “denoised” and frequency‐filtered BOLD signals of various brain regions, there are multiple techniques to compute functional connectivity and Pearson correlation may not be the optimal method. Importantly, nodal community structures may vary somewhat depending on the connectivity metric used (Mohanty, Sethares, Nair, & Prabhakaran, 2020). Nevertheless, although several alternative approaches—including wavelet and frequentist coherence, mutual information, dynamic connectivity, and fractional distance metrics—have been proposed, it is as yet unclear which methods yield the most stable and biologically plausible network configurations, and we have opted to take a more conventional approach. Third, although the statistical threshold () used to generate our fAMs was conservative and established via a rational and iterative approach, it may not represent an optimal choice, and should be more rigorously established in future studies as it determines the inclusion or exclusion of a connection. Similarly, our thresholding approach to the structural connectivity matrices (based on mean nodal degree), and the sAM produced, may be sub‐optimal. Comparisons between this method and an alternative thresholding method, based on network sparsity, indicates that the former method retained several structural connections which were not extant after thresholding with the latter, with minor differences in the consensus network (computed as HC fAM ST fAM sAM: see Supplementary Figure S14). Future studies should leverage thresholding methods, such as network‐based statistics (Zalesky, Fornito, & Bullmore, 2010) which control not only for connection‐wise false discovery rate but also for familywise error rates. Moreover, as we did not ascertain the effective connectivity—the influence of vestibular nodes on affective nodes, for example—our understanding of these interactions remains limited. This, too, represents a compelling direction for further research. Fourth, our current findings must be reevaluated in the context of larger sample sizes, which were inadequate to robustly assess the effects of age and days since injury due to rank deficiency during regression. Although the regression coefficients for these terms was <0.01, indicating a minor influence on network connectivity, these weights may potentially change with a larger sample size and the inclusion of a sufficient number of subjects with injury <4 weeks (acute PCVD) or >12 weeks (chronic PCVD) from the date of scan. Notably, however, our normalized effects and differences, which are directly comparable to a Cohen's d statistic, fall in the moderate‐to‐strong range of effect sizes even for this modest sample. Finally, it is possible that differences between control and PCVD participants with respect to affective region connectivity may reflect differences in intrinsic risk‐taking behavior rather than, or in addition to, the effects of vestibular injury and recovery, and future studies should evaluate this possibility.

CONCLUSIONS

Vestibular disorders such as PPPD and PCVD are characterized not only by deficits in orientation and static and dynamic balance but also by differential recruitment of cortical sensory and integrative areas as well as regions, such as the anterior cingulate, insula, and orbitofrontal cortex, which modulate affect. A comprehensive understanding of higher vestibular function and dysfunction therefore necessitates an integrative description across the cognitive–affective, balance, and oculomotor domains. Our current findings, which were derived from analysis of resting‐state connectivities and conditioned upon underlying anatomical structure, suggest that many of these “ancillary” areas—including anterior cingulate, insula, and prefrontal cortex—are in fact tightly incorporated with the conventional components of the vestibular cortex, including the PIVC, MSO, sensory and motor cortex, and frontal and supplementary eye fields, yielding a richer description of the vestibular neuromatrix.

CONFLICT OF INTEREST

The authors declare no potential conflicts of interest.

AUTHOR CONTRIBUTIONS

Jeremy L. Smith: Methodology, formal analysis, data curation, writing – original draft preparation and reviewing and editing, and visualization. Anna Trofimova: Investigation and writing – reviewing and editing. Vishwadeep Ahluwalia: Formal analysis, investigation, resources, data curation, and writing – reviewing and editing. Jose J. Casado Garrido: Formal analysis, writing – reviewing and editing. Julia Hurtado, Rachael Frank, and April Hodge: Investigation. Russell K. Gore and Jason W. Allen: Conceptualization, supervision, project administration, funding acquisition, and writing – reviewing and editing.

Appendix S1: Supporting information

Click here for additional data file.