Connectivity-based localization of human hypothalamic nuclei in functional images of standard voxel size

Despite their critical roles in autonomic functions, individual hypothalamic nuclei have not been extensively investigated in humans using functional magnetic resonance imaging, partly due to the difficulty in resolving individual nuclei contained in the small structure of the hypothalamus. Areal parcellation analyses enable discrimination of individual hypothalamic nuclei but require a higher spatial resolution, which necessitates long scanning time or large amounts of data to compensate for the low signal-to-noise ratio in 3T or 1.5T scanners. In this study, we present analytic procedures to estimate likely locations of individual nuclei in the standard 2-mm resolution based on our higher resolution dataset. The spatial profiles of functional connectivity with the cerebral cortex for each nucleus in the medial hypothalamus were calculated using our higher resolution dataset. Voxels in the hypothalamus in standard resolution images from the Human Connectome Project (HCP) database that predominantly shared connectivity profiles with the same nucleus were subsequently identified. Voxels representing individual nuclei, as identified with the analytic procedures, were reproducible across 20 HCP datasets of 20 subjects each. Furthermore, the identified voxels were spatially separate. These results suggest that these analytic procedures are capable of refining voxels that represent individual hypothalamic nuclei in standard resolution. Our results highlight the potential utility of these procedures in various settings such as patient studies, where lengthy scans are infeasible.

In our previous study ( Osada et al., 2017 ), we employed functional images with relatively small voxels (1.25 mm isovoxels) to reduce the voxels in the boundaries between the nuclei, enabling the discrimination of individual nuclei in the medial hypothalamus using parcellation analyses. However, scanning small-voxel images requires a substantial amount of time to compensate for the low signal-to-noise ratio of the images, and lengthy scans should be avoided in many situations including clinical studies. The present study aimed to estimate likely locations of individual nuclei in the medial hypothalamus using functional images with standard spatial resolution of 2 mm isotropic voxels. Functional connectivity between individual nuclei and the entire cerebral cortex was calculated from our highly-sampled higher resolution dataset to generate cerebrocortical correlation maps for individual nuclei. The cortical correlation maps for individual nuclei were subsequently utilized to estimate the locations of hypothalamic nuclei in the standardresolution voxels in the Human Connectome Project (HCP) database  that were connected with the cerebral cortex in a spatially similar manner. Twenty datasets of 20 subjects each were obtained from the HCP database to evaluate the reproducibility of the analytic procedures.

Subjects in the higher resolution dataset
Ten right-handed subjects (6 men and 4 women, age: 27.0 ± 7.7 years (mean ± standard deviation [SD]) ranging from 20 to 39 years) participated in the scanning of higher resolution images. Written informed consent was obtained from all subjects according to the Declaration of Helsinki. The experimental procedures were approved by the Institutional Review Board of Juntendo University School of Medicine. Data and materials are available at Dryad ( doi.org/10.5061/dryad.3j9kd51g0 ), except for raw image data because sharing of raw image data was not included in the informed consent.

MRI procedures for the higher resolution dataset
All MRI data were acquired using a 3-T MRI scanner (Siemens Skyra, Erlangen, Germany). T1-weighted structural images were obtained for anatomical reference (resolution = 0.8 × 0.8 × 0.8 mm 3 ). Higher resolution functional images of 1.25 mm isotropic voxels were obtained using multiband gradient-echo echo-planar sequences (time of repetition [TR] = 4.0 sec, echo time [TE] = 41.6 msec, flip angle = 73°, field of view [FOV] = 160 × 160 mm 2 , matrix size = 128 × 128, 120 contiguous slices, and multiband factor = 4) ( Ogawa et al., 2018 ). TR of 4 sec is more vulnerable to aliasing, and frequency higher than the filter also contains connectivity information ( Chen and Glover, 2015 ;Lewis et al., 2016 ). The small FOV (160 × 160 mm 2 ) did not always cover the whole brain along the anterior-posterior axis; the aliasing artifacts in the frontal lobe from the occipital lobe were minimal, and the aliasing overlap artifact observed in the frontal lobe was 0.1%-4.9% of the volume in the whole cerebral cortex. We acquired 100 volumes in each fMRI run at the resting state and repeated it for 10 runs in each of the 10 daily sessions, and 10,000 total volumes were collected for each of the 10 subjects. Runs of 100 images were used to inspect our subjects frequently to keep them from falling asleep. The net time taken to collect the total higher resolution images of one subject in the present study (10,000 × 4 sec = 40,000 sec) was approximately 12 times longer than that in the HCP database ( Van Essen et al., 2013 ) (4,800 × 0.72 sec = 3,456 sec). More scanning parameters are provided in Table S1.

HCP datasets
Functional images of the standard spatial resolution of 2 mm isotropic voxels were taken from the HCP database (TR = 0.72 sec, TE = 33.1 msec, flip angle = 52°, FOV = 208 × 180 mm 2 , matrix size = 104 × 90, 72 contiguous slices, multiband factor = 8) ( Van . Twenty datasets of 20 subjects each (400 subjects in total) were used to detect individual hypothalamic nuclei. Since the HCP data were acquired with the phase encoding in the left-to-right and right-to-left directions, unilateral hypothalamus (right and left, respectively) was severely distorted, whereas the hypothalamus in the other hemisphere (left and right, respectively) appeared almost intact. We therefore only used the intact side of the hypothalamus. Frame-wise displacement (FD) ( Power et al., 2012 ), a measurement of instantaneous head motion, was used for quality control of the resting-state data. Subjects were excluded from the analysis when the functional images of less than 2,000 images (out of 4,800 images in total scanned for each subject) were available after frames with FD > 0.25 mm were censored, as well as uncensored segments of data lasting fewer than 5 contiguous volumes. Sixty-nine subjects out of 400 (17.25%) were excluded from analysis.

Overview of image data analysis
Two types of resting-state functional images (1.25 mm and 2 mm) were analyzed ( Fig. 1 ). The hypothalamus was parcellated for each nucleus in the higher resolution images, and a region of interest (ROI) was defined in each subject for each hypothalamic nucleus ( Fig. 2 ). The correlation map of the cerebral cortex was generated from each hypothalamic nucleus, by calculating the correlation between the nucleus ROI and the cortical parcels distributed by Glasser et al. (2016) ( Fig. 3 ). Then, in the standard resolution images, the voxel in the hypothalamus was identified, the cortical correlation map from which was the most similar to the cortical correlation map in the higher resolution images for each hypothalamic nucleus ( Figs. 4 and 5 ).

ROIs defined for individual hypothalamic nuclei in the higher resolution
The higher resolution functional images were preprocessed for the resting-state functional connectivity, as described previously ( Ogawa et al., 2018 ;Tanaka et al., 2020 ). Briefly, images were corrected for slice timing and realigned to the first image of the first day for each subject using SPM8 ( www.fil.ion.ucl.ac.uk/spm/ ). Temporal filters (0.009 Hz < f < 0.08 Hz) were applied to the functional images using FSL ( Smith et al., 2004 ). A general linear model was used to regress out nuisance signals that correlated with head motion, whole-brain global signals, averaged ventricular signals, and averaged white matter signals. We evaluated the amount of head motion by using FD ( Power et al., 2012 ). Frames with FD > 0.2 mm were censored, as well as uncensored segments of data lasting fewer than 5 contiguous volumes. The parameter of FD = 0.25 mm for the standard resolution images was determined based on Gordon et al. (2016) , while that of FD = 0.2 mm for the higher resolution images was determined empirically keeping a sufficient amount of data, as used in our previous study ( Ogawa et al., 2018 ). The included images were 75.1 ± 12.6% (mean ± SD) of the total acquired images. The standard HCP preprocessing pipeline was not used for the higher resolution images because the small FOV (160 × 160 mm 2 ) did not always cover the whole cerebral cortex (see MRI procedures).
Parcellation analyses based on boundary mapping ( Cohen et al., 2008 ) were applied to the hypothalamus in the higher resolution images. Each voxel in the hypothalamus of each subject was used as a seed to calculate its correlations with the voxels in the gray matter of the cerebral cortex. Voxel-wise correlation coefficients in the cerebral cortex were transformed to Fisher's z. The similarity of the spatial patterns of the correlation maps was then evaluated using correlation coefficients, and similarity maps were generated. The similarity map during boundary mapping for each subject was spatially normalized to the Montreal Neurological Institute (MNI) template, using parameters (linear transformation, trilinear interpolation) between individual T1 images co-registered with functional images (affine transformation) and the standard MNI T1 image. The normalized similarity maps during boundary mapping were further deformed linearly such that the locations of the anterior commissure (AC) and the MB in the functional images, which were visually identifiable landmarks, were matched to those in the normalized T1weighted structural image ( Osada et al., 2017 ) (Fig. S1). The deformed maps were interpolated to a resolution of 1 × 1 × 1 mm 3 . After spatial smoothing (full width at half-maximum [FWHM] = 2 mm) of the Overview of image data analysis. Two types of resting-state functional images (1.25 mm and 2 mm) were analyzed. The hypothalamus was parcellated for each nucleus in the higher resolution images. The correlation map of the cerebral cortex was generated from each hypothalamic nucleus by calculating correlation between the hypothalamic ROI and the cortical parcels. Then, in the standard resolution images, the voxel in the hypothalamus was estimated, the cortical correlation map from which was the most similar to the cortical correlation map in the higher resolution images for each hypothalamic nucleus.
hypothalamus, excluding the third ventricle of each side, spatial gradients of the similarity maps during boundary mapping computed for each seed voxel. A three-dimensional watershed algorithm ( Vincent and Soille, 1991 ) was applied to the gradient maps, and the binary watershed maps were averaged across seed voxels to generate a boundary probability map. High probability indicates that the boundaries of nuclei are likely to exist and low probability indicates that the centers of nuclei are likely to exist.
Each hypothalamic nucleus in each subject was defined as the region in the hypothalamus within a 1.5 mm radius around the voxel with the smallest boundary probability, with reference to the coordinates of the hypothalamic nuclei reported in Osada et al. (2017) . Voxels in the third ventricle were excluded from the ROIs. The defined ROIs were AH (anterior nucleus of the hypothalamus), ARC (arcuate nucleus of the hypothalamus), DMH (dorsomedial nucleus of the hypothalamus), MPO (medial preoptic nucleus), PH (posterior hypothalamic nucleus), PVH (paraventricular nucleus of the hypothalamus), and VMH (ventromedial nucleus of the hypothalamus).

Correlation maps of the cerebral cortex in the higher resolution
Each hypothalamic nucleus has its specific profile of functional connectivity with the cerebral cortex. To generate such cortical correlation maps, the image data were projected to 32k fs_LR surface space using the multimodal surface matching method (MSMSulc) ( Robinson et al., 2014( Robinson et al., , 2018Glasser et al., 2016 ). We assumed that the approximate leftright symmetry in the hypothalamus, and the ROI of each nucleus in the hypothalamus consisted of voxels in the left and right sides. The time series in the ROI of each hypothalamic nucleus were extracted to calculate the correlation coefficients with the 360 cerebrocortical parcels provided in Glasser et al. (2016) . The correlation coefficients in the cortical parcels were transformed to Fisher's z. Since the small FOV did not always cover the whole cerebral cortex (see MRI procedures), parcels near the frontal or occipital pole did not always have correlation values. Out of 360 parcels in the cerebral cortex ( Glasser et al., 2016 ), 2 parcels lacked data from 5 subjects, 1 parcel lacked data from 4 subject, 4 parcels lacked data from 3 subjects, 11 parcels lacked data from 2 subjects, and 11 parcels lacked data from 1 subject. The cortical correlation maps were then averaged across subjects, and all the cortical parcels in the group-average correlation map had z values.

Correlation maps of the cerebral cortex in the standard resolution
The standard resolution functional images taken from HCP database were preprocessed for the resting-state functional connectivity, based on the HCP pipeline ( Glasser et al., 2013 ). Functional images were realigned, topup distortion corrected and spatially normalized to the MNI template. Time series data were cleaned using the ICA + FIX method ( Salimi-Khorshidi et al., 2014 ). The image data were projected to 32k fs_LR surface space using the multimodal surface matching method (MS-MAll) ( Robinson et al., 2014( Robinson et al., , 2018Glasser et al., 2016 ). Since the HCP greyordinate format contained only partial volumes of the hypothalamus, the volumetrically preprocessed images were used for the hypothalamus. Global signals were then regressed out using mean greyordinate time series.
Cortical correlation maps in the standard resolution were generated that would later be used to compare with the group-average cortical correlation map for each hypothalamic nucleus in the higher resolution. The time series was extracted from each voxel in the hypothalamus (after minimal spatial smoothing of 2 mm in FWHM), and the correlation coefficients were calculated with the 360 cerebrocortical parcels provided in Glasser et al. (2016) . In a separate analysis, we also included the subcortical structures in the correlation maps. In the HCP, the subcortical structures were defined as the nucleus accumbens, amygdala, brain stem, caudate, cerebellum, ventral diencephalon, hippocampus, pallidum, putamen, and thalamus ( Glasser et al., 2013 ). Therefore, the whole brain consisted of 360 cortical parcels and 19 subcortical regions. To avoid distorted regions of the functional images in the hypothalamus, only the intact side of the images were used for analysis, calculating the left hypothalamus-whole cortex correlation and the right hypothalamuswhole cortex correlation separately. The correlation coefficients in the cortical parcels were transformed to Fisher's z. We again assumed the approximate left-right symmetry in the hypothalamus, and the correlation maps were averaged across the left and right sides.

Similarity maps of the hypothalamic nuclei in the standard resolution
The similarity was then estimated by calculating the spatial correlation between a set of cortical correlation maps from the hypothalamic voxels in the standard resolution and a group-average cortical correlation map for a particular hypothalamic nucleus in the higher resolution. The correlation of the paired images yielded one coefficient for each hypothalamic voxel, which generated a similarity map of the hypothalamus volume for one particular nucleus in the left or right hemisphere for each subject in the standard resolution. After Fisher's z transformation, the similarity map for each nucleus was averaged across the subjects in each of the 20 datasets.

Localization of hypothalamic nuclei
The hypothalamus was delimitated manually based on the anatomical landmarks in the structural and functional images such as the AC, optic chiasm, third ventricle, globus pallidus, internal capsule, and MB, as described in previous MRI studies of the hypothalamus Schönknecht et al., 2013 ; see Table 1 in Schönknecht et al., 2013 for detailed definition). A small volume was defined for each hypothalamic nucleus to estimate the location of each nucleus based on the similarity map. The centers of the small volumes were placed based on the left/right average of the Y and Z coordinates reported in Osada et al. (2017) , and the spatial extent of the small volumes in the X axis were 2 voxels (fixed at X = 0 and ± 2, where X coordinates were defined as the center of the voxels, according to the volumetric image format of the HCP database), and 3 voxels and 3 voxels in the Y and Z axes, respectively. The voxels in the cerebrospinal fluid and the MB were excluded from the small volumes. The size of the small volumes seems empirically reasonable because, given the size of the entire hypothalamus, it is hard to imagine that the position of the nuclei moved by 4 mm. In fact, the small volumes of 5 voxels and 5 voxels in the Y and Z axes covered the anterior (supraoptic)-middle (tuberal)-posterior (mammillary) distinction of the hypothalamus ( Singh, 2011 ; Baroncini et al.,  2012 ; Schönknecht et al., 2013 ) in the Y axis and the dorsal-ventral distinction in the Z axis (Fig. S2). The voxel in the small volume that showed the highest similarity in the similarity map was regarded as the most likely location of the nucleus in the standard resolution. To evaluate the reproducibility of the estimated location of the hypothalamic nuclei, 20 datasets of the HCP database were examined.

Results
The hypothalamus in the higher resolution images was parcellated for each subject, and the centers of individual hypothalamic nuclei were estimated. Fig. 2 A shows the locations of the individual nuclei in one representative subject. The same maps were also shown in coronal slices ( Fig. 2 B). The means and SDs of the centers across the subjects were presented in Fig. 2 C. The nucleus ROIs were defined for these nuclei based on these centers for subsequent analyses.
A group-average cortical correlation map for each hypothalamic nucleus was generated ( Fig. 3 A) by calculating the correlation with the cerebral cortex in a parcel-based manner. The cortical correlation maps are provided in Supplementary data, which can be visualized using Connectome Workbench software ( www.humanconnectome.org/software/connectome-workbench ) with the standard 32k fs_LR surface file (balsa.wustl.edu). The z values of the top 20 parcels are listed for each nucleus in Table 1. We then examined similarity between the correlation maps for these nuclei by calculating the spatial correlations between these correlation maps ( Fig. 3 B). A dendrogram was generated that showed the hierarchical binary tree by defining dissimilarity between each pair of the cortical correlation maps of hypothalamic nuclei as "1 -r " ( Fig. 3 C). The dendrogram shows that the hypothalamus consisted of 2 major clusters of (1) MPO, PVH, DMH, and AH and (2) VMH, ARC, and PH. MPO and PVH showed the most similar spatial pattern of cortical correlation maps.
A similarity map was then generated by calculating spatial correlation between a set of correlation maps from the hypothalamic voxels in each dataset in the standard resolution and a correlation map for a particular hypothalamic nucleus in the higher resolution. One example of the similarity maps of the entire hypothalamus (for MPO) is shown in Fig. 4 A. Higher-value voxels in the maps of the hypothalamus indicate more likely locations of the nucleus in the standard resolution. However, the high-value voxels existed not only in likely regions but in remote regions, suggesting low specificity of the estimation. We therefore defined a small volume (2 × 3 × 3 voxels in the x, y, and z axes) to estimate the most likely location of each hypothalamic nucleus, with an estimated voxel in the higher resolution placed in the center of the small volume ( Fig. 4 B).   p < 0.001). The tests revealed significance after correction for multiple comparisons across nuclei. The similarity values in the voxels, estimated as the hypothalamic nuclei, are provided in Table S2. Fig. 5 B shows the locations of the most similar voxels that appeared, as such, most frequently in the 20 datasets.
The most likely locations of the nuclei were estimated by combining the data from the left and right sides of the hypothalamus and the cerebral cortex. We also conducted the entire analyses for the left and right sides of the hypothalamus and the cerebral cortex separately. The results of the final locations, as well as the interim results, are shown in Fig. S4 for the left side and in Fig. S5 for the right side. Although they were largely similar, the results were relatively less reproducible than those from both sides. We further conducted the entire analyses with the data in the subcortical structures included in the correlation maps. As shown in Fig. S6, the inclusion of the subcortical structures generated reproducible results (see Table S3 for Fisher's z-values of the subcortical areas). Finally, we identified the datasets which yielded insufficient similarity values in the small volumes of the hypothalamic nuclei, with the threshold of Fisher's z = 0.087 (p < 0.05). Fig. 5 C shows the estimated locations in the datasets above the threshold, whereas Fig. 5 D shows those below the threshold. The low-value datasets appeared mainly in the ventral nuclei such as VMH and ARC, and the locations of VMH and ARC in the above-threshold datasets were restricted to the upper rows of the respective ROIs.

Discussion
This study presents analytic procedures that allow us to estimate likely locations of individual nuclei in the medial hypothalamus using functional images of the standard resolution (2 mm isotropic voxels). The voxels that represent individual nuclei were reproducible across 20 datasets. The estimated voxels for the individual nuclei were spatially separable in most cases. These results suggest that the spatial profiles of functional connectivity between the individual nuclei and the cerebral cortex are capable of refining the locations of the individual nuclei, and provide the opportunity for a nucleus-level investigation of the medial hypothalamus using functional images of standard spatial resolution.
The voxels representing individual nuclei were identified with the analytic procedures based on the assumption that the nuclei can be defined by the spatial profiles of functional connectivity with the cerebral cortex. Our results showed that the spatial distribution of the identified voxels was significantly reproducible. One caveat of the analyses would be that the identified voxels ultimately may need to be clarified with histological and physiological investigations ( Williams and Elmquist, 2012 ;Sternson, 2013 ;Krashes et al., 2016 ;Lo et al., 2019 ). The present analysis was based on the results of Osada et al. (2017) where it was shown that the locations of identified nuclei were consistent with those from human histological studies ( Baroncini et al., 2012 ). An invasive histological approach can hardly be applied to living humans and very little is known to date about the physiological properties of hypothalamic nuclei in noninvasive functional imaging at the nucleus level. Future animal studies combining functional MRI and invasive approaches will help to clarify the histological and physiological properties of the identified voxels.
Global signal regression is known to induce negative correlation in the correlation maps ( Murphy et al., 2009 ;Fox et al., 2009 ;Weissenbacher et al., 2009 ;Caballero-Gaudes and Reynolds, 2017 ), and it was finally concluded that there is no correct way to preprocess the resting state data ( Murphy and Fox, 2017 ). In our previous study of higher resolution data ( Osada et al., 2017 ), the global signal regression was utilized. Our present analytic procedures estimated the similarity between correlation maps in the standard and higher resolutions. Therefore, inclusion of global signal regression in calculating the correlation maps in the standard resolution would be more appropriate in our analyses.
The spatial profiles of functional connectivity for individual nuclei were calculated with the whole cerebral cortex in this study, rather than a part of the cerebral cortex limited to the areas, such as the orbitofrontal cortex and medial prefrontal cortex, that have monosynaptic connections ( Ongür et al., 2000 ;Barbas, 2007 ;Tost et al., 2010 ;Babalian et al., 2019 ). It is well known that the functional connectivity primarily reflects monosynaptic anatomical projections ( Honey et al., 2009( Honey et al., , 2010. It is also known that the functional connectivity reflects indirect projections from various connection patterns such as common input/output ( Adachi et al., 2012 ;Bassett and Sporns, 2017 ;Reid et al., 2019 ). In a separate analysis where connectivity with subcortical structures was included in the spatial profile (Fig. S6), the final results were almost unchanged, presumably because of the relatively small number of the subcortical structures (cortex: 360, subcortical structures: 19). Thus, the functional connectivity with the whole cerebral cortex and, perhaps, subcortical structures may maximize the ability to refine the locations of individual hypothalamic nuclei.
Our results showed reproducible detection of peak voxels in the similarity maps for individual nuclei across the 20 HCP datasets, in that the peaks were clustered in the small volumes ( Fig. 5 A). The clusters in the small volumes may reflect the structural extent of the nuclei or inter-dataset variability of the locations of the nuclei. An atlas shows that the sizes of the nuclei are not confined to a 2 mm isotropic voxel ( Mai et al., 2016 ). Therefore, our analytic procedures do not have the precision to narrow them down to exactly one voxel. Importantly, when the values of the similarity between correlation maps in the standard and higher resolutions were small ( Fig. 5 ), the estimated positions of VMH and ARC should be interpreted carefully.
This study presents the results from one case of scanning condition in the standard resting-state dataset of HCP. Since the distortion of the functional images varies depending on the scanning parameters and the scanners themselves, it would not be practically sound to utilize a kind of probabilistic map of locations of the hypothalamic nuclei that can be compared with other groups scanned under different scanning conditions. The analytic procedures, instead, can be applied to functional images scanned with various scanning parameters, with a constraint of a voxel size of approximately 2 mm, which is smaller than most of hypothalamic nuclei. The present study proposes collecting sufficient resting-state functional images at approximately 60 minutes as seen typically in HCP ( Glasser et al., 2013 ;Gordon et al., 2017a ). This is done in approximately 20 normal subjects with estimating locations of the hypothalamic nuclei using analytic procedures; where our correlation maps in the higher resolution ( Fig. 3 A), available as surface files in the Supplementary data, were referred to for similarity evaluation. Then, the estimated locations can be applied to other groups of subjects, such as patient groups, scanned with the same scanner with the same parameters.