Introduction

Since its first description by Auguste Forel in 1877, the zona incerta (ZI) has remained one of the most elusive regions of the human brain. Located in the subthalamic region, the ZI is situated between the ventral thalamus, subthalamic nucleus (STN), and the red nucleus. This small gray matter region is surrounded by the H1 and H2 fields of Forel, the H field, as well as a number of important white matter pathways including the corticospinal, cerebellothalamic, and medial lemniscal tracts (Lau et al., 2020; Morel et al., 1997). Recent studies demonstrate that the ZI exhibits broad connections with the cortex, subcortical structures, brain stem, and spinal cord, resulting in theories regarding the role of the ZI as an integrative hub across a broad range of functions (Arena et al., 2024; Haber et al., 2023; Mitrofanis, 2005; Ossowska, 2020), including fear processing (Lin et al., 2023; Venkataraman et al., 2021; Zhou et al., 2021, 2018), sleep (Liu et al., 2017, 2011), attention (Chometton et al., 2021; Mitrofanis et al., 2004), and locomotion regulation (Garau et al., 2023; Kim et al., 2023; Liu et al., 2011; Sharma et al., 2018).

Based on observations in animal models, the ZI has been identified to consist of a unique molecular environment that is rich in terms of cellular and functional diversity (Chometton et al., 2017; Moriya and Kuwaki, 2020; Wu et al., 2023; Yang et al., 2022). The ZI has been determined to be predominantly composed of inhibitory GABAergic neurons with additional populations such as glutamatergic, dopaminergic, and melanin-concentrating hormone neurons, and upwards of 20 different neuron populations (Fratzl and Hofer, 2022; Mitrofanis et al., 2004; Watson et al., 2014). Based on these cytoarchitectural profiles, the ZI can be divided into 6 internal neuronal subregions in rodents: pars rostropolaris, or rostral ZI (rZI), pars dorsalis, or dorsal (dZI), pars ventralis, or ventral (vZI), pars magnocellularis, pars retropolaris, and pars caudalis (Kawana and Watanabe, 1981; Romanowski et al., 1985). The latter 3 sectors are later collectively known as the caudal ZI (cZI) (Arena et al., 2024; Nicolelis et al., 1995), whereas the dZI and vZI are commonly referred to as the central region in contemporary literature (Monosov et al., 2022; Watson et al., 2014). Furthermore, these ZI subregions exhibit differential connectivity patterns with the cortex, supporting a critical role of cortico-incertal connections for brain-wide communication across a range of functions (Barthó et al., 2002; Lin et al., 2023; Venkataraman et al., 2021; Wang et al., 2024; L. Zhang et al., 2022; Zhang and van den Pol, 2017; Zhou et al., 2021). The extent to which rodents recapitulate findings in the human ZI remains poorly understood due to substantial phylogenetic differences. Not only is the human cortex much larger in size, but it also has regions with no corresponding homologues in mice. Genetic differences, such as primate-specific neural progenitors driving cortical development, and distinct functional connectivities further highlight cross-species disparities. Therefore, direct research using primate data is essential to fully understand the ZI in humans.

Poor direct visualization of the ZI has contributed to the challenges with performing complementary analyses in humans. The cZI remains the most well-studied subregion in humans given its efficacy as a deep brain stimulation (DBS) target for Essential Tremor (ET) (Awad et al., 2020; Plaha et al., 2011; Shapson-Coe et al., 2021) and Parkinson’s Disease (PD) (Bingham and McIntyre, 2024; Blomstedt et al., 2018; Mostofi et al., 2019; Plaha et al., 2006). However, the lack of visibility using standard magnetic resonance imaging (MRI) sequences has led the surgical target to be more conventionally referred to as the “posterior subthalamic area” given it has been unclear what structure is being targeted (Blomstedt et al., 2010). Manual segmentations in stereotactic (MNI) space provide an avenue for estimating the ZI location and, when combined with diffusion MRI (dMRI), point to growing evidence of functional duality with the rZI putatively modulating obsession-compulsive behavior and the cZI modulating tremor (Saluja et al., 2024). Studies in non-human primates have provided some of the critical evidence validating these findings in the rostral and caudal ZI using dMRI with validation via tract tracing methods (Haber et al., 2023). These studies build evidence for dMRI as a viable method for characterizing cortical-incertal connections in vivo.

Building on the aforementioned work, we aim to provide a refined analysis of the topographic organization of the human ZI by mapping its cortical, i.e., cortico-incertal, connections. Our approach includes both “classical” discrete clustering via spectral clustering and the use of diffusion map embedding to identify continuous connectivity gradients. Diffusion map embedding, a method for mapping transitions in structural connectivity, provides nuanced views of connectivity differences that might resemble those observed in animal models better (Margulies et al., 2016; Royer et al., 2024). We combine both computational neuroanatomy techniques with our previous in vivo labeling of the human ZI, based on ultra-high field (7 Tesla; 7T) MRI T₁ mapping (Lau et al., 2020), and high-quality dMRI datasets from the Human Connectome Project (HCP) (Van Essen et al., 2012). Our structural connectivity workflow replicates established protocols that have previously been used to study other subcortical regions (Lambert et al., 2012). We demonstrate that our ZI structural connectivity maps (1) recapitulate known subregions from previous rodent and primate studies, (2) are replicable across field strengths and individual subjects, (3) can be linked to previously established cortical hierarchies and functional specializations, and (4) offer a proof-of-concept for personalizing neuromodulation strategies.

Results

Characterizing human cortico-incertal connectivity using in vivo tractography

We reconstructed cortico-incertal structural connectivity using in vivo probabilistic tractography based on high-resolution 7T dMRI HCP data (N=169) (Van Essen et al., 2012) and our previously validated ZI segmentation (Figure 1a-b) (Lau et al., 2020). The connections between the cortex and ZI were predominantly between the ZI and frontal lobe, representing approximately 65-75%, and secondarily connections were identified with the parietal lobe at 15-25%, with contributions identified to the temporal and occipital lobe were much more sparse (<5%) (Supplementary Figure 1). The majority of connections between the ZI and frontal lobe passed by way of the internal capsule.

Figure 1

To investigate whether a more granular topography of cortico-incertal organization was present, the initial analysis of the group-average tractogram was refined by applying diffusion map embedding and spectral clustering to the group-average structural connectivity matrix. For each subject, cortico-incertal connectivity was represented as M-by-N matrices, where M corresponds to the voxels in the ZI (N=1981 for the left hemisphere and N=1901 for the right hemisphere) and N corresponds to the cortical parcels from the HCP multimodal parcellation (HCP-MMP1.0, N=180 per hemisphere). These matrices detail the streamline counts between the ZI seed voxels and cortical brain regions (Figure 1c). Diffusion map embedding and spectral clustering of the connectivity matrix, averaged across all subjects, subsequently yielded continuous gradients of connectivity differences among ZI voxels (ranked based on their explained variance) and discrete clusters of ZI voxels with similar structural connectivity patterns, respectively (Figure 1d-f).

Mapping continuous gradients of cortico-incertal connectivity reveals a predominant rostral-caudal organization of the ZI

The diffusion map embedding results demonstrated a predominant differentiation along the rostral-caudal axis of the ZI (Figures 2a and 2b). The primary gradient (G1) extracted using diffusion map embedding distinguished between the rZI (blue) and cZI (red), while the secondary gradient (G2) identified the interface between them (red) and the rest (blue). G1 reflected a connectivity gradient from voxels highly connected to primary sensorimotor (high G1) to more anterior prefrontal cortical areas including the frontal pole (low G2, Figure 2c). G2 captured an additional gradient of ZI voxels connected to premotor and dorsolateral prefrontal areas (low G2) versus the rest of the brain (high G2). Note that additional gradients were computed (see Supplementary Figure 2a).

Figure 2

Identifying discrete ZI subregions with differential cortico-incertal connectivity patterns

In parallel with the gradient approach, we selected a spectral clustering solution of k=6, that integrates cytoarchitectonic knowledge of the number of ZI sectors identified in rodent models (i.e., rZI, dZI, vZI and cZI) (Romanowski et al., 1985) with the added benefit of providing a more granular representation of the ZI given use cases include stereotactic targeting (see final section). Other cluster solutions were also computed (k=2-8, Supplementary Figure 2b). A composite visualization of the cluster-wise tractograms, illustrating the rostral-caudal organization of cortico-incertal connections, is displayed in Supplementary Figure 3, and their top 5% of streamlines are summarized in Supplementary Table 1.

Spectral clusters 1 to 6 followed the change in G1 values, where Cluster 1 represented the lowest G1 values and Cluster 6 the highest (Figure 2b). Looking at the cluster-wise tractograms (Figure 3), the rZI involving Clusters 1-2 were demonstrated to be structurally connected with the anterior and ventral prefrontal cortex including the frontal pole (BA10), orbitofrontal cortex (BA47), as well as the pregenual and subgenual components of the anterior cingulate cortex (ACC), consistent with studies using tract tracing in macaques and dMRI in humans (Haber et al., 2023; Saluja et al., 2024). The cZI (Clusters 5-6) was demonstrated to be highly connected to the primary motor cortex (BA04; shown in green in Supplementary Figure 1b), consistent with known literature in rodent models (Yang et al., 2022) as well as human and non-human primates (Haber et al., 2023; Saluja et al., 2024).

Figure 3

Using this analysis and similar to observations in rodent and non-human primate literature (Haber et al., 2023; Yang et al., 2022), we identified a central region of the ZI (Clusters 3-4) between the rostral (Clusters 1-2) and caudal (Cluster 5-6) regions that was demonstrated to be highly connected to the intervening dorsal prefrontal cortex (BA06, BA08) with decreasing frontopolar connections moving more caudally and increasing primary motor and sensory connections.

Linking continuous and discrete connectivity findings

Having demonstrated the topographic organization of the cortico-incertal connectivity using diffusion map embedding as well as spectral clustering, how do the results using these two methods compare? To address this question, we used G1 and G2 to formalize a “2D gradient coordinate space” for the ZI, driven by the differences in cortical connectivity among individual voxels (Supplementary Figure 4a) and demonstrate the voxel-wise correspondence with the spectral clustering results. Combining G1 and G2 values (explaining 37% and 14% of the variance on average across left and right hemispheres, Supplementary Figure 2) in a 2D gradient coordinate space showed clear separation of the k=6 spectral clustering solution (Supplementary Figure 4b). This is reflected by an average Silhouette score of 0.324 (95% CI [0.311 – 0.337]) and 0.284 (95% CI [0.269 – 0.299]) for left and right hemispheres, which ranges from -1, indicating incorrect clustering, to 1, indicating highly dense clustering. Scores around zero indicate overlapping clusters. The 2D gradient coordinates were transferred to the ZI voxel and cortical surface space using a 2D colormap, and thus reduced in dimension to a single color code. This 1D color coding confirmed the dominant rostral-caudal axis among ZI voxels and their corresponding differences in connectivity from anterior prefrontal cortex to primary sensorimotor cortices (Supplementary Figures 4c and d).

Tractographic evidence of a dorsal-ventral organization of the ZI

While our current results confirmed a predominant rostral-caudal topographic organization, evidence from animal models would suggest the presence of an additional dorsal-ventral organization (Mitrofanis, 2005; Monosov et al., 2022; Watson et al., 2014). Given the potential significance of this complementary organization, we specifically investigated coronally oriented reconstructions to assess for potential differences in cortical connectivity along this axis (Supplementary Figure 5) in comparison to the in vivo T₁ maps (Lau et al., 2020), given these differences would not be easily detectable in the axial plane. Coronal cross-sections demonstrated high correspondence of G2 with the T₁-based segmentation of the internal structures of the ZI seed label (in blue and magenta), based on the Schaltenbrand histological atlas (Lau et al., 2020; Morel et al., 1997; Schaltenbrand, 1977). A full range of coronal cross-sections is displayed in Supplementary Figure 6.

Validating human cortico-incertal connectivity findings across datasets

Given that much of our analyses relied on high-resolution dMRI data obtained with optimized sequences on 7T MRI machines (Vu et al., 2015), which may not be widely accessible, it was important to assess the broader applicability of our findings. To ensure the technical validity of our results, we conducted several validation analyses. These included evaluating the replicability of our findings across images acquired with a lower magnetic field strength, assessing the reliability between test-retest sessions of the computational solutions, and confirming replicability at the individual subject level.

Replicability between MRI field strengths and test-retest reliability

Visual inspection of the group-average gradients indicated overall high replicability and reliability across field strengths (3T and 7T) and the 3T test-retest datasets, respectively. The most notable differences were between 7T and 3T (Figure 4a-b). To quantify their similarity, Procrustes disparity was calculated for each dataset pair (where lower disparity indicates higher similarity). Quiver plots in Supplementary Figure 7 visually illustrate the shift in G1 and G2 values between datasets for each voxel in both the 2D coordinate as well as volume space. Quantitative analysis confirmed the highest similarity among 3T datasets (Figure 4c) versus those compared to the 7T dataset (average disparity score of 0.065 ± 0.016 vs. 0.014 ± 0.003 across both hemispheres, p < 0.001).

Figure 4

Similar to the gradient data, the clustering results showed higher similarity among 3T datasets, as indicated by higher Dice coefficients and lower centroid distances (Figure 5). However, an average Dice score of 0.767 ± 0.042 and distance of 0.900 mm ± 0.199 mm between 7T and other datasets indicated good overlap or similarity. Notably, centroid distance results revealed a hemisphere-specific trend, with larger differences observed in the right hemisphere. Cluster-wise evaluation in Supplementary Figure 8 showed that in particular clusters 2 to 4 had the lowest replicability and reliability while clusters 5 and 6 had the highest.

Figure 5

Replicability at the individual subject level

Assessing the individual subject level replicability of the diffusion map embedding and clustering results helped to demonstrate the potential downstream utility of these techniques. High individual subject level replicability is a requisite for potential clinical applications such as DBS targeting. We found a high correlation between the group-based gradients and those from the individual subject’s data (mean Pearson’s correlation coefficient = 0.924 ± 0.052, Figure 6a). However, the primary gradient showed significantly higher values (0.963 ± 0.016) compared to the second gradient (0.886 ± 0.048, p < 0.001), consistent across both hemispheres. Compared to the between dataset comparisons (e.g., 3T vs. 7T), individual subject level cluster solutions relative to group-based clusters showed greater divergence, with centroid distances exceeding 1 mm (2.456 ± 0.246 mm) and Dice scores below 0.4 (0.228 ± 0.113, Figure 6b-c). Additionally, slightly larger differences were observed for the right hemisphere, as indicated by centroid distances, which varied among clusters more than Dice similarity.

Figure 6

Linking cortico-incertal connectivity to cortical functions and hierarchies

After confirming the reliability of the structural connectivity estimates, we focused on understanding how the observed topographic organization relates to the functional role of the ZI. To achieve this, we compared the observed cortico-incertal connectivity patterns (Figure 2c) with known cortical maps of cognitive functioning and hierarchical organization, providing a context based on extent brain maps regarding the structural and functional implications of coricto-incertal connections. Examples of these cortical properties are shown in Supplementary Figure 9.

Connectivity gradients

The cortical expression of G1 showed significant correlations with cognitive terms related to motor processing including “movement” (Pearson’s r = 0.680, p < 0.005), “coordination” (Pearson’s r = 0.565, p < 0.05), and “multisensory processing” (Pearson’s r = 0.562, p < 0.05, Figure 7a). G2 correlated strongly with terms such as “pain” (Pearson’s r = 0.479, p < 0.05), “response inhibition” (Pearson’s r = 0.444, p < 0.001), and “cognitive control” (Pearson’s r = 0.438, p < 0.001). Regarding cortical hierarchies (Figure 7b), G1 exhibited the strongest correlations with the first principal component of cognitive terms from the NeuroSynth database (“CogPC1”, Pearson’s r = 0.694, p < 0.005) (Yarkoni et al., 2011) and the magnetoencephalography theta-band activity (5–7 Hz) over the frontal regions (Pearson’s r = 0.821, p < 0.001) (Shafiei et al., 2022; Van Essen et al., 2012).

Figure 7

Spectral clustering

Strong correlations with cognitive terms and cortical hierarchies were particularly evident for clusters 3-5 (Figure 7a-b), which have the highest number of connecting streamlines (Figure 3c). Cluster 5, located near the central sulcus, is significantly linked to movement-related cognitive processing (Pearson’s r = 0.623, p < 0.005) and CogPC1 (Pearson’s r = 0.593, p < 0.005). Cluster 4, situated more anteriorly, overlaps with regions involved in working memory (Pearson’s r = 0.494, p < 0.001) with a high level of expression of serotonin (5-HT1B) receptors (Pearson’s r = 0.413, p < 0.005) (Beliveau et al., 2017). Cluster 3, located further towards the frontal pole, is associated with mood (Pearson’s r = 0.473, p < 0.005) and impulsivity (Pearson’s r = 0.473, p < 0.001), and the sensorimotor-association axis (Pearson’s r = 0.583, p < 0.005) (Sydnor et al., 2021). Clusters 1, 2 and 6, characterized by the least number of connecting streamlines (Figure 3c), were relatively weakly associated (i.e., low Spearman’s coefficient and/or within the spatial autocorrelation range) with cognitive terms or cortical hierarchies.

Clinical relevance of cortico-incertal connectivity for stereotactic targeting

Finally, to evaluate the implications of data-driven ZI structural connectivity characterization for stereotactic targeting, we assessed the overlap of the stimulation volume obtained in a patient with essential tremor (ET) treated successfully using cZI DBS, with ZI voxels in their 2D gradient coordinate space and clusters (Figure 8a-b). This patient notably had complete resolution of tremor symptoms based on the essential tremor rating assessment scale at one year follow-up post-DBS with absence of side effects, and thus could aid in defining optimal stimulation targets based on the cortico-incertal structural connectivity. For the left hemisphere, ZI voxels within the stimulation volume were characterized by mean G1 and G2 scores of 0.066 ± 0.058 and -0.130 ± 0.036, respectively, and overlapped predominantly with cluster 4 (∼64%, Figure 8c). The right hemisphere stimulation value was characterized by mean G1 and G2 scores of 0.174 ± 0.045 and -0.046 ± 0.061, respectively, and overlapped with clusters 4 (∼22%) and 5 (∼20%).

Figure 8

Discussion

Traditionally understood as “a space between structures”, our understanding of the ZI has evolved to recognize distinct subregions connected with the neocortex, which provide the anatomical basis for location-dependent functional variability (Chometton et al., 2017; Z. Li et al., 2021; Lin et al., 2023; Mitrofanis et al., 2004; Moriya and Kuwaki, 2020). In this study, we investigated the cortico-incertal connections in humans through probabilistic tractography, reporting differential connectivity patterns and their functional implications within the ZI. One key finding stemming from this approach is the identification of a topographic organization that follows a rostral-caudal axis, with evidence for a central ZI region situated between the better characterized rostral and caudal regions of the human ZI. We demonstrate that gradient and spectral clustering approaches provide complementary strategies for elucidating the connectivity of substructures within the ZI. Combined with information extracted from contextual analyses, these findings offer valuable insight regarding the structural organization of the ZI and how it relates to function and targeted therapy.

Cortico-incertal connectivity follows a topographic organization in humans

The human ZI exhibits extensive connections throughout the neocortex with predominant connectivity to the frontal lobe, followed by the parietal lobe, and relatively restricted connectivity with the temporal and occipital lobes (Supplementary Figure 1; Supplementary Table 1). Probabilistic tractography further revealed a topographic organization of the cortico-incertal connections along the rostral-caudal axis (Figure 2), validated by replicability and reliability analyses from both the 3T and 7T datasets (Figures 4 and 5). Specifically, we observed that the rZI exhibits robust connectivity with the anterior and ventral PFC as well as the rostral and subgenual ACC. Conversely, the cZI demonstrates prominent connectivity with the primary motor, premotor, and somatosensory cortices. These findings mirror evidence from previous tract tracing studies in rodents (Mitrofanis and Mikuletic, 1999; Yang et al., 2022) and non-human primates (Coudé et al., 2018; Haber et al., 2023; Ongür et al., 1998; Saluja et al., 2024), demonstrating separable connectivity patterns between rZI and cZI subregions. Specifically, rodent studies show the rZI receives inputs from the frontal lobe and the prelimbic areas, while the cZI receives inputs from the motor areas and the basal ganglia (Yang et al., 2022). Non-human primate studies revealed a similar topographic pattern, with the rZI receiving projections from the prefrontal cortex and the cZI connecting to motor areas (Haber et al., 2023; Ongür et al., 1998).

The primary gradient follows a rostral-caudal pattern of organization

A rostral-caudal cortico-incertal connectivity pattern is directly visible in the primary gradient (G1, Figure 2). Specifically, the rZI exhibits extensive streamlines anchored to the prefrontal cortex, including the anterior prefrontal cortex (BA10) and the anterior cingulate cortex, whereas the cZI is preferentially anchored with the primary sensory and motor cortices (Figure 2). These results are in line with previous tractography and axon tracing experiments (Haber et al., 2023; Saluja et al., 2024). The rostral and caudal ZI were connected to associated cortical regions by way of the anterior and posterior limbs of the internal capsule, respectively providing a structural substrate for the functional diversity observed across ZI subregions, which at the cortical level are significantly associated with coordination, sensory integration, and movement as determined via contextual analyses (Figure 7). These results are also in line with theories regarding the rostral-caudal organization of the granular PFC, from the perspective of a hierarchy of cognition that becomes increasingly more complex more rostrally, before integrating with emotional control from the agranular ventromedial and orbitofrontal cortices (Badre and D’Esposito, 2009; Koechlin and Summerfield, 2007; Levy, 2024). Previous studies in mice report that inputs from the prefrontal region to the rZI are associated with regulating escape speed (Chou et al., 2018; Hormigo et al., 2023), curiosity (Ahmadlou et al., 2021; Monosov et al., 2022), fear (Z. Li et al., 2021; Venkataraman et al., 2021; L. Zhang et al., 2022; Zhou et al., 2021, 2018), and anxiety processing (Z. Li et al., 2021; L. Zhang et al., 2022; Zhou et al., 2021). Consistent with rodent literature, recent primate studies provide evidence that the rZI acts as an integrative hub that connects cortical and subcortical targets modulating fast, survival-based responses in a dynamic environment (Haber et al., 2023). Similarly, previous findings also indicate a connection between the sensorimotor cortex and cZI, highlighting its role in modulating motor control. Selective activation of cZI neurons has been shown to alleviate motor symptoms in Parkinsonian mice (L. X. Li et al., 2021). Tract tracing experiments have further elucidated projections from the motor cortex to the ZI in non-human primates (Haber et al., 2023; May and Basso, 2018; Saluja et al., 2024). Taken together, the primary gradient reported in this study aligns with current knowledge of the rostral and caudal ZI’s involvement in diverse neurological processes, highlighting the distinct connectivity patterns of these ZI subregions to their respective cortical regions as the structural basis for their functional roles.

The secondary gradient reveals a central region in the human ZI

A second gradient (G2) was identified that is anchored on connections between the central ZI and the dorsal PFC separable from cortico-incertal connections with the rostral and caudal ZI respectively (Figure 2). To our knowledge, this is the first report of separable cortical connections to the central ZI region in humans. Supporting these results, tract tracing experiments in rodents have reported the central sectors exhibit distinct cortical connectivity patterns compared to the rostral and caudal regions (Yang et al., 2022). Detailed macaque tract tracing work has also identified clear connections between the ACC and dorsal PFC and the central ZI both dorsally and ventrally (Haber et al., 2023). G2 was found to be significantly associated with terms such as “pain”, “cognitive control”, and “response inhibition”, suggesting an integrative role of the ZI in pain processing and decision making. In line with this, optogenetic and chemogenetic experiments demonstrate that selective activation of the ZI reduces pain perception in mice (Hu et al., 2019; Li et al., 2023, 2022; Moriya and Kuwaki, 2020; Wu et al., 2023) and rats (Moon et al., 2016), whereas selective inhibition increases pain perception (Moon et al., 2016; Moriya et al., 2020; Wu et al., 2023). DBS of the ZI has been demonstrated to alleviate thermal-related pain in PD patients (Lu et al., 2021), suggesting a role for pain-related modulation.

The pattern of G2 aligns with previous findings of high resolution longitudinal T₁ mapping in the rZI (Lau et al., 2020), as voxels with higher G2 values correspond to longer T₁ relaxation time (Supplementary Figure 5). This highlights a dorsal-ventral organization that has been captured in rodent (Mitrofanis, 2005; Paxinos and Watson, 2006; Watson et al., 2014) and non-human primate models (FitzGibbon et al., 2000; Haber et al., 2023; Hardman and Ashwell, 2012; May and Basso, 2018; Mitrofanis et al., 2004). Previously, we demonstrated that T₁ mapping serves as a robust method for visualization of the ZI in vivo at 7T, separating the rostral and caudal subregions from surrounding white matter tracts (Lau et al., 2020). Building on this work, we report that probabilistic tractography permits detailed visualization of the ZI in humans, separating the rostral subregion from the surrounding fields of Forel (Lau et al., 2020), identifying a consistent organization of cortico-incertal connections.

Discrete subregions of the zona incerta using spectral clustering

While diffusion map embedding uncovered smooth transitions of cortical connections within the ZI (Margulies et al., 2016; Royer et al., 2024), mapping discrete subregions of the ZI to discrete cortical regions is a conventional, and furthermore intuitive, approach as well (Petersen et al., 2024). To explore this alternative view, we assessed the topographic organization of the ZI using spectral clustering, identifying boundaries between clusters of voxels with similar connectivity patterns. Cluster-wise findings are potentially more relatable to tract-tracing literature (Donahue et al., 2016) and may offer higher certainty compared to voxel-wise analysis, which try to resolve features at a sub-areal scale (Petersen et al., 2024).

Inspired by the cytoarchitectonic classification of the rodent ZI into six sectors (Romanowski et al., 1985), we employed a spectral clustering solution (k=6) to account for this known diversity, which furthermore provides a level of granularity for investigating location-dependent neuromodulation along the ZI (Awad et al., 2020; Bingham and McIntyre, 2024; Blomstedt et al., 2018; Mostofi et al., 2019; Plaha et al., 2011, 2006; Shepherd et al., 2024). The resulting clusters followed the same rostral-caudal topographic organization observed with diffusion map embedding (Figure 2, Figure 3). Our spectral clustering results show that rostral clusters (Clusters 1-2) are preferentially connected to cortical areas located more at the frontopolar and ventral prefrontal cortex, thus associated with cognitive and emotional processing (e.g., “facial recognition,” “mood,” “impulsivity,” “autobiographical memory”). The centrally located Clusters 3-4 connect with dorsal PFC regions, contextually linked to memory-related terms (“working memory,” “maintenance,” “rehearsal”), while caudal clusters (Clusters 5-6) connected with premotor, and M1 and S1 are associated with sensorimotor functions (“locomotion,” “coordination,” “multi-sensory”). This spectral clustering approach was found to provide evidence complementary to the gradient approach for revealing connections between premotor areas and a central ZI region, separable from the rostral ZI, which connects to the prefrontal cortex, and the cZI, which connects to the primary motor region. While the association between the rZI with cognition (Haber et al., 2023; Z. Li et al., 2021; Lin et al., 2023) and the cZI with locomotion are well established (Awad et al., 2020; L. X. Li et al., 2021), our contextual analysis reveals that the central ZI regions are implicated in prefrontally-mediated memory processing in humans. Supporting this, the ZI has been shown in mice to be essential for memory formation employing methods involving selective inhibition (Sun et al., 2024; Zhou et al., 2018).

Enhancing stereotactic targeting of the ZI through topographic mapping

Mapping the connectivity profile along the rostral-caudal axis of the ZI can help elucidate its integrative role as a neuromodulation target for a wide spectrum of neurological disorders (Awad et al., 2020; Bingham and McIntyre, 2024; Blomstedt et al., 2018; Mostofi et al., 2019; Plaha et al., 2011, 2006; Shepherd et al., 2024). The structural connectivity gradients and cluster solutions derived in this study, provide a reference framework to project DBS electrode coordinates (Figure 8). We demonstrate how these topographic maps can be employed for visualizing stimulation sites in the context of cortico-incertal connections. Compared to spectral clustering, the presented 2D gradient coordinate space provides a unique visualization of ZI organization that could be employed for refining targeting, optimizing the therapeutic window, and complementing recent developments in this area of research to identify DBS sweet spots (Hollunder et al., 2024).

Challenges of in vivo tractography

Integrating dMRI findings with validation methods such as tract tracing from the literature is crucial for assessing the accuracy of connectivity studies (Thomas et al., 2014). Our observations align well with previous tract tracing evidence from the rodent or non-human primate studies, supporting the reliability of dMRI and tractography in exploring cortico-subcortical connectivity reliably (Kai et al., 2022). However, it is important to acknowledge the spatial limitations of anatomical fidelity in these techniques (Maffei et al., 2022; Maier-Hein et al., 2017).

Identifying fine tracts, such as the dorsal-ventral organization within the ZI, is constrained by the quality of dMRI data, including spatial resolution and signal-to-noise-ratio (SNR). To address these challenges, we used the 7T HCP dataset, which was acquired with MRI protocols optimized for high spatial resolution (1.05 mm³) and SNR dMRI data (Vu et al., 2015). Furthermore, by focussing primarily on larger (cortico-incertal) white matter tracts, the inherent risk of spurious tracts is mitigated, and our validation analyses support this, showing replicability at lower magnetic field strengths and at the level of individual subjects (Figures 4-6). Note that the individual replicability metrics are sensitive to the quality of the spatial transformations between image spaces (e.g., native vs. template). Integration of highly sensitive metrics of image coregistration quality, for example, by the placement of anatomical fiducials on structural MRI scans, can be considered in future work (Lau et al., 2019; Taha et al., 2023). Alternatively, individualized ROIs can be considered to prevent the need for coregistrations between different image spaces. Moreover, the lower individual replicability observed for the spectral clustering results (i.e., relative to the gradients) are likely due to the stronger impact of single voxel differences on the respective evaluation metrics.

Tractography relies on assumptions regarding water diffusion directionality in tissues that may potentially lead to inaccuracies, particularly when dealing with complex fiber orientations as observed in subcortical regions (Axer et al., 2011; Basser and Pierpaoli, 1996). Each step in the processing workflow involves decisions that impact subsequent steps and results (Jeurissen et al., 2019; F. Zhang et al., 2022). For example, mask selection influences the number of streamlines identified; larger masks yield more streamlines, while smaller masks may miss connections. Mask accuracy is crucial, as overlap with gray matter or cerebrospinal fluid can cause anatomically implausible connections. Our ZI mask was based on extensive intra- and inter-rater reliability assessments limiting the potential for inaccuracy in labeling (Lau et al., 2020). The choice of tractography algorithm is another important consideration. Deterministic algorithms may result in more false-negatives (missed connections), while probabilistic algorithms may produce more false-positives (implausible connections). Additionally, the parameters for defining tractography termination criteria must balance sensitivity and specificity.

Concluding remarks

Using computational neuroanatomy approaches grounded in probabilistic tractography, we aimed to improve the understanding of the topographic organization of the human ZI and its connections with the cerebral cortex. Our findings highlight a distinct rostral-caudal gradient that links specific cortical areas to corresponding regions within the ZI. This detailed mapping, validated across multiple datasets and at the individual subject level, supports the involvement of the ZI across a diverse range of functions, from motor control to cognitive control and emotional regulation. The complementary evidence provided by gradient and spectral clustering approaches enhance our understanding of the structural organization of the ZI with potential applications for optimizing targeted neuromodulation strategies.

Materials and methods

Datasets

Minimally pre-processed data as part of the Human Connectome Project (HCP) young adults study (Glasser et al., 2013; Van Essen et al., 2012) were used to evaluate the structural connectivity of the ZI through probabilistic tractography. Structural connectivity analyses using the dMRI data were carried out on three separate datasets: “7T” (i.e., the reference), “3T” and “3T test-retest” (N=36; 11M/25F, aged 22-35). All subjects included in the 7T dataset (N=169, 68M/101F; aged 22-35) were included in the 3T dataset (N=169, 68M/101F; aged 22-35).

For each subject, T₁-weighted (T₁w) 3T MRI scans were acquired with a 3D MPRAGE sequence (Mugler III and Brookeman, 1990): resolution = 0.7 mm isotropic voxels; repetition time/echo time (TR/TE) = 2400 / 2.14 ms. Acquisition of dMRI data varied between field strengths. The 7T dMRI images were collected on a Siemens MAGNETOM 7T MRI system (Vu et al., 2015) with a 1.05 mm3 nominal isotropic voxel size, TR=7000 ms, TE=71.2 ms, b-values=1000, 2000 s/mm2 (64 directions per shell), FOV=210 × 210 mm2 with 15 b-value = 0 s/mm² images. For the 3T and 3T test-retest datasets, dMRI data were collected with a 1.25 mm3 nominal isotropic voxel size, TR=5520 ms, TE=89.50 ms, b-values=1000, 2000, 3000 s/mm² (90 directions per shell), FOV=210 × 180 mm² with 18 b-value = 0 s/mm² images. The 3T and 3T test-retest datasets were acquired on customized Siemens Skyra 3T MRI systems. Full acquisition details are described in the HCP1200 reference manual (https://www.humanconnectome.org/study/hcp-young-adult/document/1200-subjects-data-release).

Structural connectivity

Probabilistic tractography

Probabilistic tractography was performed to derive a structural connectivity matrix from the dMRI data. As part of the minimal preprocessing pipeline data release (Glasser et al., 2013), all subjects underwent FreeSurfer processing (v5.3.0-HCP) (Fischl, 2012). The ZI ROI mask was first resampled and transformed to the individual subjects’ minimally preprocessed volume space (0.7mm³) (Lau et al., 2020). Volumetric neocortical labels were built by mapping the HCP-MMP1.0 surface parcellation using Connectome Workbench’s ribbon-constrained label-to-volume-mapping function and FreeSurfer-derived surfaces (Glasser et al., 2016; Marcus et al., 2011). The ZI mask (thresholded at 50% and dilated with a 3 voxel radius) was used to seed tractography to the target neocortical regions using FSL’s probtrackx with 5000 streamlines per ZI ROI voxel (Behrens et al., 2007). The resulting probability maps in the ZI quantified the number of streamlines that reached each target. Generated connectivity was then transformed from subject’s native space to the MNI152NLin6Asym template space and reduced to a 2-dimensional M-by-N matrix, where M represents the voxels in the ZI ROI (1981 and 1901 for left and right hemispheres, respectively) and N is the neocortical targets (180 each hemisphere) with their corresponding number of streamlines.

Connectivity gradients

Nonlinear dimension reduction was used to transform the connectivity matrices into a low-dimensional representation (Coifman et al., 2005). ZI voxels that are characterized by similar connectivity patterns will have a value closer together in the low-dimensional space, whereas voxels with little or no similarity are farther apart. A total of hundred connectivity gradients were calculated using the GradientMaps function in the BrainSpace Python toolbox (Vos de Wael et al., 2020). These were computed in two ways: using (i) the group-averaged and (ii) the individual subject’s (for individual subject level replicability analyses, see respective section) connectivity matrices as input, while using the normalized angle kernel and the diffusion map embedding approach, and otherwise default parameter settings (Vos de Wael et al., 2018).

Spectral clustering

In parallel to the gradients, spectral clustering of the structural connectivity matrices was performed using the SpectralClustering function in the scikit-learn library for Python, using the default parameters and a chosen solution (i.e., number of clusters) ranging from k=2 to k=8 (Pedregosa et al., 2011). As for the calculation of gradients, spectral clustering calculates an affinity matrix from the connectivity matrices and performs a low-dimension embedding but adds an additional K-means clustering to assign discrete labels to each voxel. As a result, voxels belonging to the same cluster have similar connectivity patterns.

Neocortical connectivity maps

We generated neocortical connectivity maps of the spectral clustering and gradients results to visualize and evaluate the ZI connectivity in terms of its relation to the neocortex and its properties. First, gradient-weighted neocortical maps were created by multiplying each row of the ZI-neocortical connectivity matrix with the corresponding gradient value of that ZI voxel (Guell et al., 2020; Katsumi et al., 2023). Second, the cortex was associated with the ZI spectral clustering result by assigning each target neocortical region with the label of the cluster demonstrating the greatest connectivity (i.e., a “winner-takes-all” based on the maximal number of streamlines). Finally, cluster-wise neocortical connectivity maps were extracted, based on the average number of streamlines connecting the voxels within each cluster label to each neocortical parcel for contextual analysis.

Validation analyses

7T vs. 3T replicability and test-retest reliability

Connectivity gradients

To enable comparison between datasets (e.g., 7T vs. 3T), all connectivity gradients (N=100) were aligned to the 7T dataset using Procrustes shape analyses. Procrustus analysis performs scaling/dilation, rotations, and reflections of the input gradients to match the reference 7T gradients and minimize the sum of the squares of the pointwise differences (i.e., disparity) between the two input datasets, A and B:

This is particularly useful in the case of comparison of ND gradient spaces, which are characterized by arbitrary scalings. The resulting, aligned gradients were then mapped back onto the ZI voxel space to visualize and further analyze structural connectivity patterns. The same analysis was used to quantify the similarity of the ZI structural connectivity patterns among datasets based on the first two gradients (i.e., 2D gradient space) by calculating their disparity scores.

Spectral clustering

For the evaluation of the spectral clustering results, we used the Dice coefficient and centroid distance for each cluster. For the Dice coefficient, voxels of the cluster are first binarized before comparing their overlap by identifying the ratio of corresponding non-zero voxels across the two datasets to the total number of non-zero voxels across both datasets:

To calculate the centroid distance, the coordinates of the centroid for each cluster is first determined before a Euclidean distance is calculated between corresponding pairs of centroids across datasets.

Together, these two measures offer insights into the similarities of corresponding labels across the datasets. A higher Dice coefficient is better and indicative of greater volumetric overlap. The Dice coefficient is typically interpreted as follows: poor (<0.2), fair (0.2 – 0.4), moderate (0.4 – 0.6), good (0.6 – 0.8) and excellent (>0.8) overlap (Kreilkamp et al., 2019). On the other hand, a shorter centroid distance is preferable, suggesting that the extent of similarity of volumes are minimal. It can also provide insight into the extent at which the volumes differ along the spatial components.

Individual subject level replicability

In addition, considering the potential application of these gradients and spectral resulting results for stereotactic targeting at the individual level, individual subject level replicability in the 7T dataset was assessed as well. We estimated the spatial correlation between individual and group-level structural connectivity for each of the two retained gradients (Yang et al., 2020). Similarly, we also calculated Dice scores and centroid distances to further evaluate the alignment between individual subject and group-level clusters (Kai et al., 2022; Lau et al., 2019).

Contextual analyses

The gradients-based and cluster-wise neocortical connectivity maps were compared with brain-wide cognitive terms and cortical hierarchies to evaluate the functional and biological relevance of the ZI structural connectivity.

First, we used the NeuroSynth meta-analytic database (neurosynth.org) to assess how the neocortical connectivity maps of the ZI gradients and cluster-specific connectomes relate to cognitive terms (Yarkoni et al., 2011). NeuroSynth aggregates data from over 14,000 published functional MRI studies to derive probabilistic measures of the association between brain voxels and cognitive processes based on frequently occurring keywords such as “pain” and “attention”. Similar to the approach in Hansen (2023) based on the Cognitive Atlas, a public ontology that provides an extensive list of neurocognitive processes (Poldrack et al., 2011), we utilized 124 cognitive and behavioral terms, covering broad categories (“attention”, “emotion”), specific cognitive processes (“visual attention”, “episodic memory”), behaviors (“eating”, “sleep”), and emotional states (“fear”, “anxiety”).

Second, we used the neuromaps toolbox to assess whether the ZI neocortical connectivity maps are linked to molecular, microstructural, electrophysiological, developmental, and other functional properties (N=73) of the neocortex (Markello et al., 2022). See “Statistical analysis” for more details.

Clinical relevance

Finally, as a proof-of-concept, we reconstructed the stimulation volumes from a single DBS case targeting the ZI to characterize the structural connectivity patterns of the voxels located within the stimulation volumes in both hemispheres.

Data acquisition

Pre-operative T₁-weighted (T₁w) and T₂-weighted (T₂w) 7T MRI scans were acquired. T₁w imaging was performed using the 3D MP2RAGE pulse sequence (Marques et al., 2010): resolution = 0.8 mm isotropic voxels; repetition time/echo time (TR/TE) = 6000/2.69 ms; inversion times = 800/2700 ms; flip angles = 4/5°. T₂w-imaging was performed using a 2D SPACE pulse sequence with parameters: resolution = 0.6 mm isotropic voxels and TR/TE = 2000/131 ms. During preprocessing, images were corrected for gradient non-linearity distortions. MP2RAGE data were additionally corrected for B ⁺-inhomogeneities using a separately acquired Sa2RAGE B ⁺ map (1.9×1.9×2.8 mm voxels; TR/TE = 2400/0.81 ms; TI1/TI2 = 45/1800 ms; flip angles = 4/11°) (Eggenschwiler et al., 2012; Haast et al., 2018; Marques et al., 2010).

Stimulation volume reconstruction

We employed Lead-DBS 3.0 (lead-dbs.org) (Horn and Kühn, 2015) for electrode localization and reconstruction in a common space (see Figure 8a). In brief, a linear co-registration was performed between the pre-operative, preprocessed T₁w MRI and post-operative CT. Subsequently, the 7T T₁w and T₂w images were non-linearly registered to the ICBM 2009b Nonlinear Asymmetric (i.e., MNI2009bAsym) (Fonov et al., 2009) template space. All registrations were performed in accordance with previously validated presets (Avants et al., 2008) of the Advanced Normalization Tools (ANTs; stnava.github.io/ANTs/) (Avants et al., 2011) as implemented in Lead-DBS. The SimBio algorithm (Oostenveld et al., 2011; Wolters et al., 2006) was used to estimate the stimulated volume based on stimulation parameters at least 1-month after the DBS device had been turned on. The essential tremor rating assessment scale scores were time-locked with DBS stimulation parameters.

Statistical analyses

For each cortical map derived from the NeuroSynth and neuromaps databases (see “Contextual analyses” section), we parcellated it using the HCP-MMP1.0 atlas and computed its Pearson’s spatial correlation with the ZI neocortical connectivity maps, accounting for spatial autocorrelations using N=10,000 spin tests. During each spin test, the parcel coordinates were randomly rotated, and original parcels were reassigned the value of the closest rotated parcel according to the Hungarian algorithm (Kuhn, 1955). Additionally statistical evaluations were performed using the appropriate tests implemented in the pingouin v0.5.3 Python package.

Acknowledgements

Data were provided in part by the Human Connectome Project, WU-Minn Consortium (Principal Investigators: David Van Essen and Kamil Ugurbil; 1U54MH091657) funded by the 16 NIH Institutes and Centers that support the NIH Blueprint for Neuroscience Research; and by the McDonnell Center for Systems Neuroscience at Washington University. Author RAMH was supported by a H2020 Marie Skłodowska Curie Actions Postdoctoral Fellowship (101061988). Author ARK was supported by the Canada Research Chairs program #950-231964, NSERC Discovery Grants RGPIN-2015-06639 and RGPIN-2023-05558, Canadian Institutes for Health Research Project grant #366062, Canada Foundation for Innovation (CFI) John R. Evans Leaders Fund project #37427, the Canada First Research Excellence Fund, and Brain Canada. Author JCL was supported by an NSERC Discovery Grant RGPIN-2023-05562 and research start-up funding from the Department of Clinical Neurological Sciences at Western University.

Data and code availability

The non-clinical data used in this study are available as part of the publicly available Human Connectome Project S1200 release (humanconnectome.org/study/hcp-young-adult). Anonymized stimulated volumes and clinical scores are available at github.com/ataha24/zona-clusters. The raw clinical data are available upon request. Analysis was performed using a combination of publicly available toolboxes: Connectome Workbench v2.0 (Marcus et al., 2011), BrainSpace v0.1.10 (Vos de Wael et al., 2020), NeuroSynth v0.3.8 (Yarkoni et al., 2011), neuromaps v0.0.3 (Markello et al., 2022), BrainStat v0.3.2 (Larivière et al., 2023), and custom Python code available at the code repository. A full list of installed Python packages is listed in the repository as well.

Additional files

Supplementary Figures 1-9

Supplementary Table 1