Towards HCP-Style Macaque Connectomes: 24-Channel 3T Multi-Array Coil, MRI Sequences and Preprocessing

Macaque monkeys are an important model species for understanding cortical organization of primates, yet tools and methods for noninvasive image acquisition (e.g. MRI RF coils and pulse sequence protocols) and image data preprocessing have lagged behind those developed for humans. To resolve the structural and functional characteristics of the relatively thin macaque cortex, high spatial, temporal, and angular resolutions are required while maintaining high signal-to-noise ratio to ensure good image quality. To address these challenges, we developed a macaque 24-channel receive coil for 3-T MRI with parallel imaging capabilities. This coil enabled adaptation of the Human Connectome Project (HCP) image acquisition protocols to the macaque brain. We also adapted HCP preprocessing methods optimized for the macaque brain, including spatial minimal preprocessing of structural, functional MRI (fMRI), and diffusion MRI (dMRI). The coil provided high signal-to-noise ratio and high efficiency in data acquisition, allowing four- and five-fold acceleration for dMRI and fMRI, respectively. Automated parcellation of cortex, reconstruction of cortical surface, removal of artefacts and nuisance signals in fMRI, and distortion correction of dMRI performed well, and the overall quality of basic neurobiological measures was comparable with those for the HCP. The resulting HCP-style in vivo macaque MRI data show considerable promise for analyzing cortical architecture and functional and structural connectivity using advanced methods that have previously only been available for humans. Highlights ➢ 24-channel 3T MR receive coil designed for the smaller macaque brain. ➢ In vivo macaque imaging protocols adapted according to guidelines from the HCP. ➢ Parallel imaging yields five- and four-fold acceleration in fMRI and dMRI sampling. ➢ HCP’s minimal preprocessing and denoising pipelines adapted for macaques. ➢ The multi-modal MRI data show considerable promise for HCP-style analyses.


Abstract: 43
Macaque monkeys are an important model species for understanding cortical organization of 44 primates, yet tools and methods for noninvasive image acquisition (e.g. MRI RF coils and pulse 45 sequence protocols) and image data preprocessing have lagged behind those developed for humans. 46 To resolve the structural and functional characteristics of the relatively thin macaque cortex, high 47 spatial, temporal, and angular resolutions are required while maintaining high signal-to-noise ratio 48 to ensure good image quality. To address these challenges, we developed a macaque 24-channel 49 receive coil for 3-T MRI with parallel imaging capabilities. This coil enabled adaptation of the Human 50 Connectome Project (HCP) image acquisition protocols to the macaque brain. We also adapted HCP 51 preprocessing methods optimized for the macaque brain, including spatial minimal preprocessing of 52 structural, functional MRI (fMRI), and diffusion MRI (dMRI). The coil provided high signal-to-noise 53 ratio and high efficiency in data acquisition, allowing four-and five-fold acceleration for dMRI and 54 fMRI, respectively. Automated parcellation of cortex, reconstruction of cortical surface, removal of 55 artefacts and nuisance signals in fMRI, and distortion correction of dMRI performed well, and the 56 overall quality of basic neurobiological measures was comparable with those for the HCP. The 57 resulting HCP-style in vivo macaque MRI data show considerable promise for analyzing cortical 58 architecture and functional and structural connectivity using advanced methods that have previously 59 only been available for humans. 60 61 Introduction 69 Old World monkeys are an important neuroscientific model for understanding primate 70 neuroanatomy (Brodmann K., 1905;Felleman and Van Essen, 1991;Van Essen et al., 2001). Macaque 71 monkeys have provided insights about neurovascular coupling (Goense and Logothetis, 2008), neural 72 wiring (Markov et al., 2014) and the evolution of the human brain's functional connectome 73 (Passingham, 2009;Wang et al., 2012). However, macaques are separated from humans by 25 74 million years of evolution, and are known to have substantial brain differences despite being 75 members of the same primate order. Recent imaging studies have revealed substantial 76 neuroanatomical differences between macaques and humans, for example in language connectivity 77 or proportion of cortex devoted to lightly myelinated association areas (Donahue et al., 2018;78 Glasser et al., 2014;Rilling et al., 2008). At the level of cortical areas, high confidence homologies 79 (i.e., a common evolutionary origin) have only been firmly established for a modest number of early 80 sensory and motor areas (Van Essen and Dierker, 2007) but are more challenging to delineate for 81 higher cognitive regions such as prefrontal cortex (Mars et al., 2018b(Mars et al., , 2018a. Improvements to in in 82 vivo neuroimaging acquisition and preprocessing may help address several outstanding questions: 83 what is the optimal interspecies registration between macaque and human cerebral cortices? What 84 are the optimal methods for non-invasively estimating functional and structural connectivity as 85 assessed by comparison with gold standard invasive tracers in macaques? What brain networks are 86 shared and which ones are different between macaques and humans? 87 88 Recently, the Human Connectome Project (HCP) developed an improved, integrated approach to 89 brain imaging acquisition, analysis, and data sharing (Glasser et al., 2016b). The overall goal of this 90 approach is to increase the sensitivity and precision with which brain imaging studies are conducted 91 in the hope that this will yield results that are more neurobiologically interpretable and more 92 accurately comparable across individuals and studies. The HCP-style approach has seven core tenets 93 (Glasser et al., 2016b): 1) Acquire as much high-quality data from as many subjects as possible. 2) 94 Acquire data with maximum feasible resolution in space and time 3) Preserve high data quality 95 throughout preprocessing by removing physical distortions, subject movement within and between 96 scans, image intensity inhomogeneities, and artefacts and nuisance signals without blurring the data 97 or altering the neural signals ( structures . 5) Align brain areas across subjects, not cortical folds (Robinson et 101 al., 2018(Robinson et 101 al., , 2014. 6) Use a data-driven structurally and functionally sensible parcellation, ideally 102 multi-modal MRI measures such as those acquired by HCP. Achieving comparable results in 126 macaques requires not only higher resolution and SNR but also low geometric distortion and signal 127 intensity inhomogeneity, and requires optimized hardware, sequences, and post-processing 128 techniques. 129

130
In this study, we designed and built a 24 channel receive coil with a geometry optimized for parallel 131 imaging of anesthetized macaque monkeys at 3T. Capitalizing on the accelerated imaging capabilities 132 of the coil, we adapted HCP-style data acquisition protocols for structural MRI , 133 fMRI (Smith et al., 2013) and diffusion MRI (dMRI) Uğurbil et al., 2013) to 134 the small size of the macaque brain, as well as the HCP-style minimal spatial preprocessing and 135 denoising pipelines (Andersson and Sotiropoulos, 2016;Glasser et al., , 2016aGlasser et al., , 2013Salimi-136 Khorshidi et al., 2014. We generate accurate white and pial cortical surfaces, subcortical 137 segmentations, myelin maps, and cortical thickness maps from structural MRI, surface aligned fMRI 138 dense timeseries that have spatial artefacts and nuisance signals removed, resting state functional 139 networks, and diffusion-based fiber orientation estimates, example tractography connections, and 140 cortical neurite orientation and dispersion imaging (NODDI) (Zhang et al., 2012 animal's MRI data was used to delineate the contour of the head surface and imported into the 3D 163 digital design software where the inner surface of the coil was designed to closely fit the surface of 164 the head (Fig. 1A). Next, 16 pentagonal and 8 hexagonal elements were configured over the surface 165 ( Fig. 1B), resembling a soccer-ball coil design (Wiggins et al., 2006). These elements were arranged in 166 three quasi-horizontal arrays to maximize parallel encoding power of multiband EPI sequences for 167 animals placed in the supine position and axial slices. The inner body of the device was constructed 168 using a 3D printer (M200, Zortrax, Olsztyn, Poland) (Fig. 1C), and the coil elements were arranged 169 over its external surface. Initially the coil elements were wired using a thin copper foil-plate (width 5 170 mm); however, because the plate elements markedly interfered with B1 transmission (data not 171 shown), the coil elements were rewired using thin coaxial copper cables (Fig. 1D,  overlap each other to reduce coupling between nearest-neighbor coils (Roemer et al., 1990), and 175 those in the caudal-posterior part were designed to have relatively larger diameter (35% larger in 176 maximum diameter) to increase sensitivity to distant brain regions (e.g. cerebellum) while reducing 177 sensitivity to closer regions (e.g. occipital cortex). The two elements placed over the eyes were also 178 relatively large in diameter to allow video recording of eyes and eyelids for monitoring depth of 179 anesthesia. In addition, capacitors were arranged vertically against the surface of the coil frame to 180 reduce interaction with B1-transmission (Fig. 1D). Fig. 1E shows the circuits, which followed a 181 standard design (Wiggins et al., 2006) consisting of diode detuning trap, cable trap and bias T 182 connected to low input-impedance preamplifiers (Siemens Healthcare, Erlangen, Germany). The 183 completed coil is shown in Fig. 1F. 184

Coil Evaluation 186
Coil elements were assessed for the ratio of loaded to unloaded quality factor Q, nearest-neighbor 187 coupling, and active detuning. Element coupling was also estimated with gradient off-line noise 188 correlation measurements. Two phantoms (NaCl 0.9%, gadolinium 0.1 mM) were designed and 189 prepared using a 3D printer: one to closely match the inner-surface of the coil (Fig. 1G) used for B1 190 quality evaluation and the other to match to a typical macaque brain size used for geometry-191 dependent noise amplification. B1-transmission was assessed with a vendor provided flip-angle 192 sequence. B1-receive field was estimated using a gradient-echo sequence and by calculating the 193 signal ratio between 24-channel and body receive coils. Finally, geometry-dependent noise 194 amplification due to parallel imaging was evaluated using gradient-echo imaging and GeneRalized 195 Autocalibrating Partial Parallel Acquisition (GRAPPA) (Griswold et al., 2002)  FOV=128x128x112mm, matrix=256×256, slice per slab=224, coronal orientation, readout direction I 232 to S, phase oversampling=15%, TR=3200 ms, TE=562 ms, bandwidth=723 Hz/pixel, no fat 233 suppression, GRAPPA=2, turbo factor=314, echo train length=1201 ms and pre-scan normalization). 234 The total acquisition time for structural scans was 22 min (17 min for T1w and 5 min for T2w). 235 236 Functional Acquisition Protocol 237 To reduce susceptibility induced geometric distortions and signal loss, the data was acquired in LR 238 and RL directions. Functional scans were acquired using gradient-echo EPI (FOV=95x95 mm, 239 matrix=76×76, 1.25 mm isotropic, interleaved slice order, and number of slices=50 covering the 240 whole brain). 241 242 An empirical estimate of the effect of multiband slice acceleration factor on fMRI tSNR was obtained 243 by a procedure similar to that used by the HCP (Smith et al., 2013). Briefly, simultaneous slice 244 excitation enables a multiband factor fold reduction in the TR and subsequent incomplete T1-245 recovery leads to a reduction in the optimal (Ernst) flip angle and thus in tSNR. However, as more 246 data volumes can be acquired in a matched time window, a more relevant estimate for the data 247 quality can be calculated by multiplying the tSNR with a square root of acquired data timepoints. following imaging parameters were applied: FOV=90 mm, matrix=100×100, 0.9 mm isotropic 285 resolution, number of slices=60, interleaved slice acquisition, multiband factor=2, GRAPPA= 2, 286 TR=3400 ms, flip-angle=90, TE=73 ms, 6/8 phase partial Fourier, echo spacing=1.12 ms, 287 bandwidth=1086 Hz/pixel, pre-scan normalization on and fat suppression using gradient reversal 288 technique (Gomori et al., 1988). Total acquisition time was 30 min, during which frequency drift was 289 small (≈0.5 Hz/min). By applying slice and in-plane accelerations (2×2), the acquisition time was 290 reduced by more than 3-fold than without acceleration. However, the shortest possible TR was not 291 used, in order to preserve SNR (to allow near-complete longitudinal magnetization recovery). 292

Animal experiments 294
Macaque monkeys (mean 5380 g, range 3030-8850 g) were initially sedated with intramuscular 295 injection of dexmedetomidine (4.5 µg/kg) and ketamine (6 mg/kg). A catheter was inserted into the 296 caudal artery for blood-gas sampling, and tracheal intubation was performed for steady controlled 297 ventilation using an anesthetic ventilator (Cato, Drager, Germany). End-tidal carbon dioxide was 298 monitored and used to adjust ventilation rate (0.2 to 0.3 Hz) and end-tidal volume. After the animal 299 was fixed in an animal holder, anesthesia was maintained using intravenous dexmedetomidine (4.5 300 µg/kg/hr) and 0.6 % isoflurane via a calibrated vaporizer with a mixture of air 0.75 L/min and O2 0. Preprocessing began with the PreFreeSurfer pipeline, in which structural T1w and T2w images were 317 registered into an anterior-posterior commissural (AC-PC) alignment using a rigid body 318 transformation, non-brain structures were removed, T2w and T1w images were aligned using 319 boundary based registration (Greve and Fischl, 2009), and corrected for signal intensity 320 inhomogeneity using B1-bias field estimate. Next, data was transformed into a standard "Yerkes19" 321 macaque atlas (Donahue et al., 2018(Donahue et al., , 2016 by 12-parameter affine and nonlinear volume 322 registration using FLIRT and FNIRT FSL tools (Jenkinson et al., 2002). 323 324 Then, the FreeSurferNHP pipeline reconstructed the cortical surfaces using FreeSurfer v5.3.0-HCP 325 (Fischl, 2012). This process included conversion of data in AC-PC space to a 'fake' space with 1-mm 326 isotropic resolution in volume with a matrix of 256 in all directions, intensity correction, 327 segmentation of the brain into cortex and subcortical structures, reconstruction of the white and 328 pial surfaces and estimation of cortical folding maps and thickness. The intensity correction was 329 performed using FMRIB's Automated Segmentation Tool (FAST) (Zhang et al., 2001) followed by 330 scaling the whole brain intensity by a species-specific factor (=80). This process significantly 331 improved white and grey contrast particularly in the anterior temporal lobe as well as white surface 332 estimation, an effect that may be associated with the so-called 'anterior temporal lobe problem' in 333 pediatric brains, potentially due to less myelination in these white matter areas. We also improved 334 the subcortical parcellation training dataset for the macaque brain, and trained for 21 subcortical 335 structures: brainstem plus bilateral accumbens, amygdala, caudate, claustrum (which is not a part of 336 the default structures for human FreeSurfer), cerebellum, diencephalon, hippocampus, pallidum, 337 putamen, and thalamus (Fischl et al., 2002). The training dataset for brain mask extraction was also 338 created. After parcellating the cortical and subcortical structures with these training datasets using 339 the T1w image, the claustrum was treated as putamen, so that subsequent white surface estimation 340 accurately estimates the white surface beneath the insular cortex, as shown in the Results. The pial 341 surface was estimated using the T2w image to help exclude dura and blood vessels, similar to the 342 HCP pipeline . We modified this procedure by applying an optimized value of 343 maximal cortical thickness (=10mm in 'fake' space, 5mm in real space like the FreeSurfer default). 344 The surface and volume data in 'fake' space was transformed back into the native AC-PC space, and 345 cortical thickness was recalculated in the animals' real physical space. were mapped onto the cortical surface using an algorithm weighted towards the cortical mid-457 thickness (Fukutomi et al., 2018). 458

Results 459
Coil performance 460 Coil bench tests showed that the unloaded/loaded Q ratio of the individual coil elements were 461 approximately 215/75=2.9. This relatively low Q-ratio results from the small degree of loading and 462 small electromagnetic flux due to the small diameter of the coil elements. Decoupling between 463 adjacent elements was less than -20 dB indicating low mutual inductance between the elements. 464 This produced noise correlation coefficients averaging 0.084 (interquartile range 0.02 and 0.126) 465 with a maximum of 0.395 (see correlation matrix Fig. 2A). High noise correlation was largely 466 constrained to the nearest neighbor elements (see Fig. 1B for element geometry, see also 467 Supplementary Fig. S1 for coil channel-specific noise correlation maps). 468

469
The inverse g-factor map, a measure of coil element separation, illustrates geometry dependent 470 signal intensity variation due to parallel image reconstruction used for dMRI (Fig. 2C). A reduction 471 factor of two yields an average inverse g-ratio slightly larger than unity (1/g=1.03 ± 0.07; values 472 reported throughout text as mean ± s.d. unless otherwise specified), indicating a small noise 473 cancellation attributable to low element noise correlation and parallel image reconstruction. 474 However, larger reduction factors (R=3 and 4) yield substantial degradation of signal intensity 475 depending on geometry (Fig. 2C, D), suggesting that a maximum GRAPPA of 2 is practical for this coil. 476 (D) The boxplot shows 1/g-factor with respect to reduction factor. While geometric distortions are small with acceleration factor of 2 (1/g=1.03±0.07), further reduction yields large signal degradations. Geometric distortions were evaluated using a phantom whose contour was matched to the average macaque brain.

477
Macaque Data Quality Evaluation 478 Structural bias-field corrected T1w and T2w weighted images acquired at 500-µm resolution are 479 shown in Fig. 3A and B for an exemplar single subject. Note the good SNR and contrast of the white 480 matter to grey matter (and to CSF) throughout the brain. 481 482 Flip-angle maps indicate that the transmission was slightly higher in subcortical regions compared to 483 cortical structures (Fig. 3C), as expected. However, the surface map (Fig. 3D) Table S1). However, a relatively low cortical tSNR in 497 lateral occipito-temporal cortex was notable in macaque data (Fig. 3H), which is mainly attributable 498 to a large B0 dephasing effect (Fig. 3F). 499 FreeSurfer automated segmentation of cortical and subcortical structures using our NHPHCP 503 structural pipeline was reliable across the subjects (Supplementary Fig. S4B), and benefited from 504 additional signal intensity normalization ( Supplementary Fig. S4A, see also Supplementary Fig. S2A  505 and S2B for B1-transmit and receive fields, respectively). Inspection of pial and white matter surface 506 contours indicates that the automatic segmentation generally followed the contrast boundaries of 507 the T1w image ( Supplementary Fig. S4C) and the T2w image ( Supplementary Fig. S4D) appropriately, 508 including in challenging thin heavily myelinated regions such as early visual and somatosensory 509 cortex. The subcortical structures including claustrum, pallidum, putamen, were automatically and 510 accurately segmented by the improved subcortical atlas (Supplementary Fig. S4B). The newly added 511 intensity normalization improved the problematic estimation of the white matter surface in the 512 anterior temporal lobe ( Supplementary Fig. S5A right), which was not achieved using the default 513 intensity bias field correction ( Supplementary Fig. S5A left). The claustrum parcellation strategy also 514 improved the white matter surface just beneath the insular cortex ( Supplementary Fig. S5B, right), 515 which often resulted in 'claustrum invagination' of the white surface by the default FreeSurfer 516 ( Supplementary Fig. S5B left). The claustrum parcellation also improved myelin contrast in the 517 anterior insular area (see next paragraph). To estimate the optimum multiband factor in fMRI, we determined the relationship between tSNR x 549 sqrt(timepoints) and multiband acceleration factor (Fig. 5A) and found that tSNR x sqrt (timepoints) 550 increases up to a factor of 5, then decreases. This pattern was more evident after the data was 551 processed using ICA-based artefact removal algorithm FIX, which yielded approximately 25% 552 improvement in tSNR. In the cortical ribbon, denoised tSNR x sqrt (timepoints) is clearly highest at 553 MBF=5 (Fig. 5B). 554

555
The resting-state fMRI runs were analyzed using multi-run sICA + FIX. The resulting sICA components 556 (a total of number of components: 124 ± 29 for each animal, N=30) were manually classified as noise 557 (on average 100 ± 23 components per animal) or signal (24 ± 9 components per animal). The manual 558 classification worked well to train FIX, and the classification accuracy achieved reasonably high 559 performance ( Using RestingStateStats in HCP Pipeline Marcus et al., 2013), the variance in 568 macaque resting-state fMRI runs was divided into six categories. Fig. 6 shows their relative 569 contributions to the total signal variance (38,400 ± 13,000, N=20, see also Table S2). Relative 570 variance estimations in descending order were unstructured noise (70.0 ± 4.8%), high-pass filtered 571 noise (15.3 ± 4.5 %), structured noise (i.e. artefacts and nuisance signals, 6.0 ± 1.5%), (neural) BOLD 572 fluctuations (4.1 ± 2.3%), motion (2.9 ± 1.3%), and FIX-denoised global signal timeseries (1.0 ± 0.7%). 573 In comparison to HCP, unstructured noise accounted for a slightly larger portion in macaque (Fig. 6), 574 which mainly originates from subcortical structures (see Supplementary Fig. S6 for spatial 575 distribution of the variance categories). Furthermore, the relative BOLD contribution was smaller in 576 macaque (4.1%) in comparison to HCP (7.7 ± 2.6%). Taken together, the contrast-to-noise ratio 577 (CNR), defined as ratio between BOLD and unstructured signal, was smaller in macaque (0.21 ± 0.07) 578 than in HCP (0.37 ± 0.08), which may be due to reduced BOLD signals in the anesthetized state (see 579 section Resting-state fMRI in Discussion). 580  Fig. S7A) and after preprocessing (Supplementary Fig. S7B) 585 demonstrated that preprocessing reduced structured artefacts. The mean global timeseries (MGT) 586 also demonstrate that FIX reduced the global signal variance, which in humans is primarily related to 587 respiration after movement artefacts have been removed by sICA+FIX. MGT power spectrum 588 ( Supplementary Fig. S7C) revealed distinct peaks within the ventilation frequency range (0.25 to 0.30 589 Hz). Preprocessing effectively attenuated ventilation artefacts, but only partially attenuated the low 590 frequency, more likely neural, fluctuations (<0.1 Hz). Across subjects, the MGT variance was 2,230 ± 591 1,530 prior to preprocessing and 170 ± 110 after preprocessing ( Supplementary Fig. S7D, N=20). 592 There appears to be relatively less global physiological noise in the macaque data relative to the 593 human data (Glasser et  less than 0.2 Hz, which are typical of RSNs. A similar functional connectivity pattern was found using 605 a single greyordinate seed placed over the area 7A (Fig. 7B). Both the RSN signal components (a total 606 of 32 signals) and the dense functional connectome can be interactively viewed in Connectome 607 Workbench after downloading data from the BALSA database (https://balsa.wustl.edu/3ggwG). 608 Overall, these results demonstrate that our experimental setup enables robust functional 609 connectivity detection and analysis. Following the HCP paradigm, we used reversed left-right phase-encoding directions in dMRI 614 acquisition to reduce TE, TR and distortion and to increase SNR and angular CNR. An example of 615 image distortion and correction (axial and coronal views) is shown in Supplementary Fig. S8. Image 616 distortions are large near regions with large B0 inhomogeneity (i.e. temporal lobe, see Fig. 3E, F). 617 Nonetheless, distortion correction was accurate, albeit with some signal drop-out and degraded SNR 618 in these regions. Mean motion absolute displacement during 30-min acquisition was 0.36 ± 0.07 mm 619 (N=10), ensuring little interaction between head motion, eddy-currents and changes in static 620 magnetic field. In contrast to HCP at 3T (Uğurbil et al., 2013), we used simultaneous MB and GRAPPA 621 acceleration to reduce distortions. Inspection of temporal stability of the dMRI acquisition did not 622 reveal pronounced structural artefacts around the ventricles and basal slices ( Supplementary Fig.  623 S9), thus indicating that simultaneous MB and GRAPPA accelerations did not substantially interact 624 with physiological noise (Uğurbil et al., 2013). The dMRI quality assurance measures were similar 625 between this study and the HCP (Fig. 8). Average SNR (whole brain) was 11.6 ± 1.4 in macaque 626 (N=10) and 9.4 ± 0.9 in the HCP (N=10) (Fig. 8A). Exemplar subject data are compared in 627 Supplementary Fig. S10 and Supplementary Table S3. The CNR slightly increased towards higher b-628 values and was similar across the studies (Fig. 8B). In white matter, three crossing fibers voxels 629 (selected by thresholding at 0.05 of third fiber's volume fraction) were detected in 59% ± 7% and 630 57% ± 4% of voxels in macaque and the HCP, respectively (Fig 9D). Finally, the dispersion 631 uncertainties of 1 st , 2 nd and 3 rd fiber orientations these voxels exhibited were also similar across the 632 studies ( Fig 9E). 633 634 Figure 8. Comparison of dMRI quality measures between macaque and the HCP (blue and red bars, respectively; N=10).
Plots show whole brain SNR (A) and CNR across b-values 1000, 2000 and 3000 (B), as well as three-crossing fiber ratio (C) and dispersion uncertainties (in degree) of 1 st , 2 nd and 3 rd fiber orientations in the white matter voxels (D). Overall, the quality measures were comparable across the studies. 635 Figure 9 shows M132 parcellated cortical maps of MD (Fig. 9A), FA (Fig. 9B), NDI (Fig. 9C) and ODI 636 ( Fig. 9D) (N=6). The MD is low in the primary motor (F1) and premotor cortices (such as F2, F4, F5), 637 and primary sensory cortices including somatosensory (areas 3, 1, 2), visual (V1) and auditory 638 cortices including core, as well as intraparietal sulcus area (Fig. 9A), whereas the NDI is high in all of 639 these areas. MD and NDI were strongly anti-correlated (R=-0.75, p<0.001). The ODI was high in the 640 periphery of the V1, somatosensory area 1, auditory cortices including core (Fig. 9D) and 641 intermediate in MT and other higher visual areas. The FA was higher in the frontal and anterior 642 temporal cortices and strongly anti-correlated with ODI (R=-0.86, p<0.001). These results are 643 comparable with those observed in the HCP (Fukutomi et al., 2018). The structural connectivity 644 patterns extracted from diffusion tractography (DT) were also parcellated and explored with respect 645 to the published quantitative retrograde tracer data (Fig. 9E, F)  Here, we have presented an adaptation of the HCP's approach to multimodal MRI acquisition, 656 preprocessing, and analysis to the macaque, using the combination of a custom-made 24-channel 657 receive-coil, high-resolution parallel imaging, and the HCP-NHP preprocessing and analysis pipelines. 658 This approach yields robust estimates of cortical thickness, myelin content, and functional and 659 diffusion measures. Importantly, since the presented protocols used share similar strengths to the 660 HCP image acquisition, and the data is stored in a common geometrical framework system ('CIFTI 661 greyordinates'), we anticipate that it will facilitate direct multi-modal comparisons with an 662 Our multichannel receive coil, fabricated to closely fit a large macaque head (Fig. 1A) will allow 672 routine imaging of macaque monkeys of different species with a range of lateral muscles and head 673 sizes. The close proximity of the coil to the head allows high SNR in the brain with further SNR gains 674 in the cortex produced by the small size of the elements (Fig. 1) (Janssens et al., 2013;Wiggins et al., 675 2006). This design allowed acquisition of both T1w and T2w structural whole-brain image acquisition 676 with a 0.5mm isotropic resolution in 22 minutes (Fig. 3a, b). In conjunction with homogeneous RF 677 transmission (Fig. 3C, D), these two features enabled automatic and robust subcortical 678 segmentations and reconstructions of pial and white matter surfaces ( Supplementary Fig. S4). 679 680 Twenty-four receive elements were arranged so as to optimize efficiency of spatial encoding capability 681 in the axial slice direction (Fig. 1B, D). This geometrical arrangement yields a relatively small noise 682 correlation coefficient (0.084), which is smaller than in previous macaque multi-channel coil designs 683 such as 0.12 in a 24-channel (Gilbert et al., 2016)  To accurately map BOLD signals onto the cortical sheet, the image resolution (1.25 mm isotropic) 693 was matched with 5 th percentile of cortical thickness (Fig. 4N, O, P) to reduce the partial volume 694 effects from white matter and CSF signals , following the HCP data acquisition 695 strategy at 3T (resolution 2 mm, the 5 th percentile of human cortical thickness ). 696 The reduction from an isotropic volume of 2 mm to 1.25 mm, however, incurs a 4-fold SNR penalty. 697 Nonetheless, tSNR of fMRI in macaque (Fig. 3G, H) is superior to that in the HCP acquired with 698 comparable imaging parameters ( Supplementary Fig. S3, Table S1). This tSNR gain may be primarily 699 attributed to the close proximity to the animal and small diameter of the receive coil elements, with 700 an additional gain from relatively small bandwidth. This illustrates the power of parallel imaging to 701 overcome a physical size difference of a factor of twelve (macaque and human brain volumes are 702 approximately 100 cm 3 and 1200 cm 3 , respectively). 703 704 While informative, tSNR is not an explicit measure of fMRI sensitivity to blood flow changes induced 705 by neural activity. It is well known that variation of fMRI signal is a mixture of nuisance (e.g. motion 706 and respiration) and neural BOLD components. To obtain insight into the content of our fMRI signals, 707 we categorized different signal sources and found that neural BOLD signal explains approximately 708 4.1% of the total fMRI variance (in data grand mean scaled to 10,000; corresponding to 773 ± 438 in 709 absolute variance) in anesthetized macaque resting-state (Fig. 6). In HCP fMRI data (awake-state), 710 neural BOLD signal explains approximately 7.7% of total variance (corresponding to 4158 ± 1594 in 711 absolute variance Marcus et al., 2013)). Because the image acquisition 712 protocols and image qualities are similar across the studies ( Supplementary Fig. S3), we speculate 713 that the lower BOLD neural signal in our macaque data may be due to, 1) attenuated thalamo-714 cortical and cortico-cortical synchronization in the anesthetized state, and/or 2) a ceiling effect of 715 signals due to relatively high blood flow, oxygen extraction rate, and saturation in anesthetized 716 macaque brain (Kudomi et al., 2005). This issue may be overcome with widely used contrast agents 717 (i.e. MION) and cerebral blood volume weighted fMRI (Mandeville et al., 1998)  we demonstrated that FIX is also very successful reducing such artefacts (6.0% of total variance, Fig.  725 6) with over 98% classification accuracy (threshold at 20, Table 1) in the macaque resting-state fMRI. 726 The relative global mean variance and its reduction in macaque (1.5% before cleanup and 1.0% after 727 cleanup) is smaller in comparison to the HCP (3.2% before cleanup and 2.2% after cleanup) (Glasser 728 et al., 2018). This smaller global signal variance in anesthetized macaques can be attributed to more 729 stable global blood flow because respirations and pCO2 were regulated by mechanical ventilation 730 (Birn et al., 2006). The majority of the signal variance, however, is unstructured noise (>60%), in 731 particular at subcortical regions that are distant from the coil elements ( Supplementary Fig. S6), 732 which can be effectively reduced using parcellation and/or Wishart filtering ( replicate several of these RSNs. Taken together, from the data quality perspective, the 24-channel 740 coil yields macaque rfMRI data that can be accurately and sensitively mapped onto cortical sheet 741 and is comparable in quality with the HCP rfMRI data, whereas from the physiology perspective, we 742 must be cautious when making inferences because of the potential effects of anesthesia on both 743 neural activity and neurovascular coupling. We will explore this topic in future work on a specialized 744 coil for awake monkey imaging. 745 746 While scaling the fMRI resolution with respect to the cortical thickness is a minimum requirement to 747 accurately localize BOLD signal within the cortical sheet, another important factor is the size of 748 functional imaging voxels relative to the area of the cortical surface for identifying sharp gradient 749 ridges in FC (Glasser et al., 2016a). We found that macaque cortical grey matter surface area is 750 ≈10,100 mm 2 per hemisphere, which is close to previous estimates of 11,900 mm 2 (Chaplin et al., for the presented dMRI protocol (0.73 mm 3 / 23.000 mm 3 ≈ 3 x 10 -5 ) approximately matches 2.5 mm 775 isotropic resolution in the human white matter (16 mm 3 / 500.000 mm 3 ≈ 3 x 10 -5 ) but is an order of 776 magnitude larger than in the HCP (1.95 mm 3 / 500.000 mm 3 ≈ 4 x 10 -6 ), although a more precise 777 comparison would require investigations on features such as radii of curvature, tract and blade 778 thickness. Smaller voxel size could aid in distinguishing challenging fiber pathways, however, under 779 our experimental conditions further reduction was impractical due to gradient power and SNR 780 limitations. 781

782
To mitigate this limitation, our strategy was to acquire data with exceptionally high angular 783 resolution (500 directions) capitalizing on two-by-two acceleration (out-of-plane MB and in-plane 784 GRAPPA) enabled by the multichannel array coil. The effect of this strategy was shown in the 785 comparable sensitivity to 3 rd crossing fibers between species (Fig. 8), despite the resolution 786 limitation in macaque. A recent ex vivo macaque study used high-quality, high-field magnetic field 787 (4.7T), long data acquisition (≈27 h) postmortem and gadolinium enhanced diffusion scans to 788 demonstrate a relatively good correspondence between probabilistic tractography and quantitative 789 retrograde tracer (R=0.55-0.60) (Donahue et al., 2016). Here, we replicated a part of those results 790 (Fig. 8E, F), thus, augmenting the findings of Donahue and colleagues to in vivo applications that are 791 within practical time limitations (≈30 min). Taken together, these results suggest that the 'HCP'-style 792 dMRI data acquisition protocols are well positioned to produce quantitative tractography measures 793 that are neuroanatomically meaningful. 794

795
The high spatial resolution with respect to cortical thickness enabled us to carry out cortical surface 796 mapping of neurite properties and to provide preliminary evidence for nonuniformity in the 797 composition and distribution of neurites in macaque cerebral cortex (Fig. 9C, D). Neurite properties 798 are considered important because the density of neurites constitute basic building units (axons and 799 dendrites) of neuronal networks, while ODI provides an indicator of the heterogeneity of neurite 800 fiber orientations, a ratio between tangential and radial fibers (Fukutomi et al., 2018). We found that 801 NDI was highest in V1 and higher than average in other visual representation areas (V2, V3, V4, and 802 MT), somatosensory (1, 2, 3 and A1), motor (M1) and granular prefrontal (Fig. 9C), cortical 803 distributions resembled those of myelin contrast (Fig. 4L). ODI was high in early somatosensory, 804 auditory and visual cortices (Fig. 9D). Together, these results are in good agreement with the HCP 805 data (Fukutomi et al., 2018). Essen and . Comparison between transitions in multimodal neuroimaging contrasts, 812 such as MT myelination (Fig. 4C) and functional connectivity (Fig. 8A, E), are particularly suggestive 813 of brain area boundaries (Glasser et al., 2016a). Therefore, the approach to data collection and 814 analysis presented here provides macaque data that may aid in multi-modal parcellation of the 815 macaque and generation of structural and functional connectomes (Glasser et al., 2016b), though a 816 robust delineation of cortical areas into functionally distinct areas will assuredly require analysis of a 817 more extensive dataset. 818 819 These HCP-style macaque data also provide an attractive substrate for multi-modal registration 820 across species-in particular, macaques and humans. Just as myelin maps and resting state networks 821 are used to register across human subjects (Robinson et al., 2018(Robinson et al., , 2014, they could be used to 822 register the cerebral cortex between group averages of humans and macaques. This would allow 823 direct comparisons between human and macaque structural and functional connectivity (Mars et al., 824 2018b). That said, we expect that the gains for cross-individual registration with areal features in the 825 macaque will be less than those in humans simply because folding patterns and the relationships 826 between folds and areas are less variable in macaques than they are in humans. Additionally, this 827 approach to macaque imaging acquisition and analysis can be used to form structural and functional 828 connectomes in the macaque for comparison with invasively measured tracer datasets (Donahue et 829 al., 2016;Glasser et al., 2016b). Such validation analyses will help to determine the optimal methods 830 for forming structural and functional connectomes in human studies , where a 831 direct comparison with a gold standard is not available. Future work will also explore cross-species 832 comparisons between macaques and marmosets imaged using specialized hardware and an HCP-833 style approach. These acquisition and analysis methods can also be applied to study disease models 834 in primate species where controlled and invasive methods can be used to investigate causality and 835 plasticity of structural and functional connectomes and their importance in shaping primate 836 behavior. 837 838

Conclusions 839
A 24-channel phased-array coil for 3T was constructed and optimized for in vivo parallel imaging of 840 macaque monkey brain. The coil provided high SNR whole-brain coverage and allowed parallel 841 imaging with high speed acquisition by a five-fold and four-fold increase in functional and diffusion 842 MRI, respectively. The data acquisition strategy in combination with the HCP-NHP minimal 843 preprocessing pipelines enabled robust mapping of structural and functional properties onto surface 844 of the cortex. The presented protocols can be acquired within a single imaging session and represent 845 compelling advance in identifying multi-modal cortical topology and structural and functional 846 connectomes in the macaque. Overall, this study demonstrates that MRI studies in animals may 847 benefit from adapting the methodologies introduced by the HCP.