A Putative Multiple-Demand System in the Macaque Brain

In humans, cognitively demanding tasks of many types recruit common frontoparietal brain areas. Pervasive activation of this “multiple-demand” (MD) network suggests a core function in supporting goal-oriented behavior. A similar network might therefore be predicted in nonhuman primates that readily perform similar tasks after training. However, an MD network in nonhuman primates has not been described. Single-cell recordings from macaque frontal and parietal cortex show some similar properties to human MD fMRI responses (e.g., adaptive coding of task-relevant information). Invasive recordings, however, come from limited prespecified locations, so they do not delineate a macaque homolog of the MD system and their positioning could benefit from knowledge of where MD foci lie. Challenges of scanning behaving animals mean that few macaque fMRI studies specifically contrast levels of cognitive demand, so we sought to identify a macaque counterpart to the human MD system using fMRI connectivity in 35 rhesus macaques. Putative macaque MD regions, mapped from frontoparietal MD regions defined in humans, were found to be functionally connected under anesthesia. To further refine these regions, an iterative process was used to maximize their connectivity cross-validated across animals. Finally, whole-brain connectivity analyses identified voxels that were robustly connected to MD regions, revealing seven clusters across frontoparietal and insular cortex comparable to human MD regions and one unexpected cluster in the lateral fissure. The proposed macaque MD regions can be used to guide future electrophysiological investigation of MD neural coding and in task-based fMRI to test predictions of similar functional properties to human MD cortex. SIGNIFICANCE STATEMENT In humans, a frontoparietal “multiple-demand” (MD) brain network is recruited during a wide range of cognitively demanding tasks. Because this suggests a fundamental function, one might expect a similar network to exist in nonhuman primates, but this remains controversial. Here, we sought to identify a macaque counterpart to the human MD system using fMRI connectivity. Putative macaque MD regions were functionally connected under anesthesia and were further refined by iterative optimization. The result is a network including lateral frontal, dorsomedial frontal, and insular and inferior parietal regions closely similar to the human counterpart. The proposed macaque MD regions can be useful in guiding electrophysiological recordings or in task-based fMRI to test predictions of similar functional properties to human MD cortex.

Pervasive human MD activation during demanding tasks suggests a core role in supporting goal-oriented behavior . One might therefore expect a similar network to exist in nonhuman primates that readily perform similar tasks after training. Many human-macaque homologies in brain function and connectivity profiles have been reported across frontoparietal association cortex Denys et al., 2004;Margulies et al., 2009;Mars et al., 2011;Hutchison et al., 2012a;Sallet et al., 2013;Miranda-Dominguez et al., 2014;Hutchison et al., 2015;Neubert et al., 2015). However, a macaque counterpart of the MD network has not yet been described. Recent studies suggest that some human frontoparietal networks may not have spatial correspondents in the monkey whether based on responses to similar task events Patel et al., 2015) or comparison of restingstate networks (Mars et al., 2011;.Therefore, the existence and structure of macaque MD regions remain open questions. Although single-cell recordings from macaque frontal and parietal cortex show some similarities to human fMRI MD responses, including adaptive coding of various task components (Watanabe, 1990;Sakagami and Niki, 1994;Toth and Assad, 2002;Freedman and Assad, 2006;Roy et al., 2010;Mante et al., 2013;Stokes et al., 2013), they necessarily come from a few prespecified locations, so could greatly benefit from knowledge of where MD foci are likely to lie.
The MD system is typically defined in humans by contrasting fMRI responses across pairs of closely matched task conditions that differ only in their level of some specific cognitive demand. In any given task, this contrast will be confounded with idiosyncratic aspects of the task or stimuli, so the generality of the response is demonstrated by meta-analyses or conjunctions that show similar effects of cognitive demand across widely varying contexts, the specific demands of which might include response conflict, perceptual difficulty, working memory load, etc. (Duncan and Owen, 2000;Fedorenko et al., 2013). Unfortunately, practical constraints make this challenging in the monkey, but an alternative is suggested by coactivation of human MD regions during rest (Fox et al., 2005;Dosenbach et al., 2007;Seeley et al., 2007;Power et al., 2011). We therefore sought to identify a macaque counterpart to the human MD system using fMRI connectivity while monkeys were lightly anesthetized, based initially on standard anatomical comparative mapping. We show that a first approximation of macaque MD regions is indeed functionally connected during anesthesia. We then refined these foci using an iterative approach to find the voxels local to each region that are most connected to the rest of the network. Finally, we tested the connectivity of these optimized regions on data withheld from the optimization step (using leave-one-out cross-validation) and describe their brain-wide connectivity.

Materials and Methods
Animals. Thirty-five (28 male) healthy Rhesus macaques (Macaca mulatta) were scanned under light anesthesia using a published protocol (Mars et al., 2011;Sallet et al., 2011). Animals weighed 4.4 -13.3 kg (mean, 8.2) and were 3.3-7.9 years of age (mean, 5.3). They were housed socially in same-sex pairs or small groups, with a 12 h light/dark cycle (lights on 0700 -1900 h) and had scheduled access to water between 12 and 16 h on testing days and at all times otherwise. All experimental procedures were performed in compliance with the United Kingdom Animals (Scientific Procedures) Act of 1986. A Home Office (UK) Project License, obtained after review by the University of Oxford Animal Care and Ethical Review Committee, licensed all procedures. The housing and husbandry were in compliance with the guidelines of the European Directive (2010/63/EU) for the care and use of laboratory animals. fMRI data from subsets of these animals have contributed to previous publications (Mars et al., 2011;Sallet et al., 2011;Mars et al., 2012;Sallet et al., 2013;Neubert et al., 2015).
Data acquisition. Data were acquired on a 3 T MRI scanner with a full-sized horizontal bore using a four-channel phased array coil (H. Kolster; Windmiller Kolster Scientific) and a local transmission coil. The anesthetized monkeys were scanned in sphinx position in an MRIcompatible stereotactic frame (Crist Instrument). A total of 872-1600 (mean, 1539) volumes of echoplanar imaging were acquired per animal, with 36 axial slices, in-plane resolution 2 ϫ 2 mm, slice thickness 2 mm, no slice gap, repetition time ϭ 2 s, and echo time ϭ 19 ms. T1-weighted structurals were acquired in the same sessions at 0.5 mm resolution using a magnetization-prepared rapid gradient-echo sequence (128 slices, resolution 0.5 ϫ 0.5 ϫ 0.5 mm, slice thickness 2 mm, no slice gap, repetition time ϭ 2.5 s, and echo time ϭ 4.01 ms).
Monkeys were initially sedated 1.5-2 h before commencing the resting-state fMRI scans, using an intramuscular injection of ketamine (10 mg/kg) combined with either xylazine (0.125-0.25 mg/kg) or midazolam (0.1 mg/kg) and buprenorphine (0.01 mg/kg). Once sedated, monkeys were given atropine (0.05 mg/kg, i.v.) to reduce secretions, a nonsteroidal anti-inflammatory drug (meloxicam, 0.2 mg/ kg, i.v.) for analgesia, and an H2 receptor antagonist (ranitidine, 0.05 mg/kg, i.v.) to protect against gastric ulceration as a side effect of nonsteroidal antiinflammatory drug treatment. Local anesthetic (5% lidocaine/prilocaine cream and 2.5% bupivacaine subcutaneous) was injected around the ears to block peripheral nerve stimulation at least 15 min before positioning the macaque's head in the stereotaxic frame. An intravenous cannula was used for delivery of fluids (warmed sterile saline drip, 5 ml/h/kg) throughout the scan. Monkeys were intubated, ventilated with intermittent positive pressure to ensure a constant respiration rate during the functional scan, and placed on isoflurane (n ϭ 33) or sevoflurane (n ϭ 2) anesthesia. The lowest possible concentration of isoflurane/sevoflurane gas (mixed with medical oxygen) was used to maintain anesthesia, in accordance with veterinary instructions. The inspired isoflurane concentration was in the range of 0.8 -2.0% (mean, 1.53%). The expired isoflurane concentration was in the range of 0.8 -1.8% (mean, 1.41%). Under these conditions, resting-state networks are reliably observed Mars et al., 2011;Sallet et al., 2013;Hutchison et al., 2014). Two animals received sevoflurane instead of isoflurane at inspired concentrations of 2.1% and 2.5% and expired concentrations of 2.0% and 2.4%. Normal body temperature was maintained throughout the scan using blankets, bubble wrap, and heated wheat bags. Respiration rate, inspired and expired CO 2 , and inspired and expired isoflurane/sevoflurane concentration were monitored and recorded using VitalMonitor software (Vetronic Services). Heart rate, blood pressure, and core temperature were also monitored throughout the scan and clinical checks assessed muscle relaxation before scanning.
Preprocessing. Analysis used MATLAB (The MathWorks), SPM8 (Statistical Parametric Mapping; www.fil.ion.ucl.ac.uk/spm), FSL (FMRIB Software Library; http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/), Caret (Computerized Anatomical Reconstruction Toolkit; Van Essen, 2012; http://brainvis.wustl.edu), and automatic analysis software (Cusack et al., 2014;www.automaticanalysis.org). The structural volume was spatially normalized to the 112RM-SL template (112 rhesus macaque template in the space of the atlas of Saleem and Logothetis, 2007) and segmented into gray matter, white matter, and CSF tissue classes (McLaren et al., 2009). The first six functional volumes were discarded to allow the signal to reach equilibrium. The remainder were used to estimate movement parameters and then registered to the structural, warped to template space using the previously calculated transformation, and spatially smoothed within each maximum-   between MD ROIs after warping to macaque and a similar number of control regions; dark and light circles indicate correlations significantly greater than or less than zero, respectively ( p Ͻ 0.05, FDR corrected). Right, Warped ROIs illustrated on the F99 cortical surface. c, As in b, but after expanding ROIs to local anatomical boundaries. d, As in b and c, but after iteratively optimizing ROIs to maximize connectivity within expanded regions (tested on left-out data); ROIs are colored where they overlap in all cross-validation folds, whereas black borders show locations present in some cross-validation folds. Dark and light dots in b and c indicate connections with strength that has significantly increased or decreased, respectively, relative to the preceding ROI definition (warped Ϫ expanded Ϫ optimized; p Ͻ 0.05, FDR corrected). e, Networkhood of MD ROIs, measured as [mean connectivity between themselves] minus [mean connectivity between MD and control ROIs]. The control ROIs are either the six selected regions illustrated in the correlation matrices of b-d (red) or all atlas regions outside of the expanded anatomical masks (blue). Error bars are 95% confidence intervals for the difference from the test values marked by correspondingly colored dashed lines. Each increase in networkhood is highly significant ( p Ͻ 0.001).
likelihood tissue compartment using a Gaussian kernel with 3 mm full-width-half-maximum. A gray-matter mask was defined per subject, as voxels with gray-matter probability Ͼ0.9. More conservative white-matter and CSF tissue masks were defined per subject as voxels with tissue probability that exceeded 0.9 in that animal and in 90% of other animals. A further single mask was drawn by hand around the posterior superior sagittal sinus in template space.
Overview of analysis strategy. Identification of putative macaque multiple demand regions proceeded in four main stages, indicated by vertical arrows in Figure 1. First, established human frontoparietal MD regions were mapped to the macaque brain using a previously published transformation. Second, these were expanded to cover local anatomically defined regions. Third, the originally warped regions were iteratively optimized to find the best-connected voxels within each expanded anatomical region. Finally, whole-brain connectivity maps from the optimized regions were combined, thresholded, and parcellated. We describe these steps in more detail below.
Human MD ROIs. Human MD ROIs were based on data from Fedorenko et al. (2013). Group-level t-statistics from contrasts that isolated cognitive demand were averaged across seven tasks. To create a symmetrical volume, data from left and right hemispheres were averaged and then projected back to both hemispheres. To define ROIs around peak locations, the volume was thresholded at t Ͼ 1.5 (Fig. 1a, left); if a cluster contained multiple peaks that also survived t Ͼ 2.7, then it was subdivided, assigning each voxel to the nearest subregion at the higher threshold. Clusters outside of the frontoparietal cortex were discarded to leave seven bilateral ROIs (Fig. 1a, right) around the presupplementary motor area and anterior cingulate (preSMA/ACC), anterior insula, intraparietal sulcus (IPS), superior and inferior precentral sulcus (sPreCS, iPreCS), and middle and anterior middle frontal gyrus (MFG). Discarded clusters included two in the occipital lobe that were suspected to reflect domainspecific responses to the visual stimuli used (Fedorenko et al., 2013) and which are not associated with fluid intelligence deficits when damaged, unlike the frontoparietal foci (Woolgar et al., 2010). Subcortical clusters were also discarded for the additional practical reason that the humanto-macaque deformation and the macaque anatomical atlas used (see below) are defined on the cortical surface. Importantly, the final analysis tested connectivity from the seed ROIs to every voxel of the brain volume. Therefore, if any nodes of the network had been missed in the initial ROI selection, then they should show up in this analysis).
Warping MD ROIs from human to macaque. Human MD ROIs were mapped from Montreal Neurological Institute (MNI152) volumetric space to the human Population-Average Landmark and Surface-based (PALS) atlas and then warped to the macaque F99 surface using Caret (Van Essen and Dierker, 2007). The deformation was based on 23 landmark homologies elsewhere in the brain (Denys et al., 2004;. See Miranda-Dominguez et al. (2014) for further validation of this deformation procedure. Regions were then converted to a volume in F99 space using Caret, and warped to 112RM-SL space at 0.5 mm isotropic resolution using SPM8. Regions were dilated by up to 3 mm in any direction within a gray-matter mask derived from the 112RM-SL tissue probability maps. If a region became fragmented, then the largest cluster was retained. After warping to the macaque, the two most anterior regions were small and closely abutting, so were combined into a single region (MFG; Fig. 1b, right). Note that we will continue to name the ROIs by their source in the human brain, although this may not correspond to macaque anatomy.
Macaque anatomical template regions. Macaque anatomical regions came from the LV-FOA-PHT composite cytoarchitectonic parcellation . The name of the parcellation derives from the initials of the authors of the parcellations on which it was based (Ferry et al., 2000;Lewis and Van Essen, 2000;Paxinos et al., 2000). Ferry et al. used the parcellation from Carmichael and Price (1994). Regions on the right hemisphere cortical surface were first converted to a volume in F99 space using Caret and warped to 112RM-SL space at 0.5 mm isotropic resolution using SPM8. Clusters of six voxels or fewer were removed and regions were dilated by up to six voxels in any direction into empty zones of a gray-matter mask derived from the 112RM-SL tissue probability maps. This ensured comprehensive labeling of gray-matter voxels by expanding from the midsurface ribbon into which the surface had been interpolated and filling in the few small regions that remained unlabeled in the composite parcellation. Finally, regions were reflected to the left hemisphere.
For comparison with MD ROIs, a similar number of control ROIs were selected from the LV-FOA-PHT set to be maximally distant from both the MD regions and from each other. These were areas 10o, 14c, 29a, 36, the floor of the superior temporal area, and the ventral posterior area. Note that no set of control regions could match all characteristics of the MD ROIs and this was not our aim. The purpose of the control ROIs was merely to show that the optimization of MD ROIs (see below) tended to increase their connectivity with each other specifically, rather than their global connectivity with any brain region. We therefore also repeated this analysis using all LV-FOA-PHT regions outside of MD seed regions as a "control." Because the final analysis uses whole-brain connectivity from the optimized ROIs, it does not depend on the particular choice of control ROIs.
Connectivity analysis. Physiological confound covariates were constructed by extracting principal components of the time series from white-matter and CSF tissue masks (as many components as required to explain 99% of variance, up to a maximum of 6; mean 1.3 for whitematter and 3.1 for CSF; Behzadi et al., 2007). A further vascular component was defined as the mean time course within the superior sagittal sinus mask. The first temporal derivatives of these time courses were also included. A motion confound covariate was defined as the time course of average displacement over the expected brain volume (approximated as a sphere of radius 40 mm; Jenkinson, 1999). The first temporal derivative of this vector was also included, as well as the element-wise squares of both vectors. Discrete cosine transform covariates were added to implement a temporal band-pass filter (0.0025-0.05 Hz) as commonly used . After projecting covariates from each gray-matter time series, functional connectivity was estimated as bivariate correlations between the mean gray-matter time series from each pair of ROIs and between these ROI time series and each voxel in the brain. Group statistics on correlation values used one-sample t tests after Fisher r-to-z transformation. p-values were thresholded to control the false discovery rate (FDR) at Ͻ0.05 (Genovese et al., 2002).
Iterative optimization. To refine MD ROIs, we used an iterative approach to find the voxels in the vicinity of each region that were most connected to the rest of the network. First, local anatomical masks were defined around each candidate ROI within which the precise location of the ROI would be optimized (Fig. 1c, right). These expanded masks included any LV-FOA-PHT anatomical regions overlapping an ROI (such that the volume of intersection with the ROI was Ͼ1% of the volume of their union).
To avoid overfitting and ensure generalization across the N ϭ 35 animals, a leave-one-animal-out cross-validation scheme was adopted. The following steps were thus performed on each set of N Ϫ 1 animals. At each iteration step, each ROI was redefined as all voxels within its local anatomical mask for which the mean correlation with other MD ROIs was significantly higher than a local baseline, defined as the mean correlation between other MD ROIs and the voxels immediately bordering the mask. (Using a local baseline per region allowed the optimization to select peaks of connectivity that were especially strong within each region even if they were weaker than peaks in other regions that had globally higher correlations, e.g., due to higher signal-to-noise or closer proximity to other MD ROIs.) Significance of higher-than-baseline correlation was assessed using one-sample t tests across the N Ϫ 1 animals. After adjusting each ROI in this way, the process was iterated until the ROIs no longer changed (Fig. 1d, right). Connectivity of the iteratively optimized ROIs was then tested using independent data from the left-out animal. Finally, group statistics were calculated across these 35 sets of crossvalidated connectivity measures.
The labels of the six MD ROIs were preserved across this optimization process. However, this partition of the optimized ROIs was not necessarily expected to be anatomically meaningful. The labels were initially retained for the purposes of comparing connectivity across initial and optimized ROIs and for visualizing the result of the optimization process. In subsequent whole-brain connectivity analyses (see below), the initial labels were abandoned and new spatial clusters were derived from the set of voxels delivered by the optimization step.
Combining, thresholding, and partitioning whole-brain connectivity maps. Original MD ROIs had been constrained to frontoparietal cortex and their subsequent optimization had been constrained within local anatomical zones. We ultimately wished, however, to identify loci anywhere in the brain that were functionally connected to this network. Whole-brain connectivity maps from each optimized ROI were therefore created for each animal using seed ROIs optimized on independent data from the other animals. Connectivity maps were then averaged across animals and positive correlations tested via one-sample t tests.
Statistical maps from the six seed ROIs were combined using the Simes method for dependent data and thresholded at an FDR of 0.025 (Heller et al., 2007) to produce a map showing the number of MD ROIs with which each voxel was significantly positively correlated. This conjunction map was largely symmetrical, so data from left and right hemispheres were averaged and then projected back to both hemispheres. The resultant volume was then thresholded and partitioned (similarly to the original human MD conjunction of cognitive demand activations) to identify locally contiguous clusters of voxels that were significantly connected to at least five of the six MD ROIs. If a connected cluster contained multiple local maxima of mean connectivity, then it was subdivided, assigning each voxel to the nearest local maximum. Clusters were discarded if they did not contain any voxels significantly connected to all six MD ROIs (two very small clusters of two and eight voxels), leaving eight bilateral clusters.
Comparison with previously reported anatomical projections. We matched all regions in the LV-FOA-PHT parcellation to brain sites in the CoCoMac (collation of connectivity data on the macaque brain) database (http://cocomac.g-node.org; Stephan et al., 2001;Bakker et al., 2012;. We then queried CoCoMac 2.0 using default options plus the S.O.R.T algorithm (Bezgin et al., 2008) to generate a list of axonal projections between these regions and to convert this to a connectivity matrix indicating projections as present (in either direction and regardless of strength), absent, or unknown. We identified the seven atlas regions that had maximum overlap with our final MD functional connectivity clusters, excluding primary auditory cortex, and assessed the proportion of known projections involving these nodes, ignoring self-connections and links for which connectivity was unknown. We also examined projections to a set of posterior parietal regions that overlapped substantially with the parietal MD regions initially warped from humans, but not with the final optimized macaque clusters. These were LV00a-LIPv, LV00a-VIPl, LV00a-VIPm, LV00a-MIP, LV00a-PO, and LV00a-5V.

Anatomically constrained optimization increases connectivity of MD ROIs
Deformation of human ROIs using landmark homologies elsewhere in the brain (Denys et al., 2004) provided a first approximation of potential macaque MD regions. It is not expected to be perfect, especially in association cortex, which tends to be far from the landmarks and has undergone the largest evolutionary expansion (Hutchison et al., 2012a). To refine these ROIs, we therefore iteratively adjusted the size/location of each MD ROI to maximize connectivity to other MD ROIs. Optimized ROIs were constrained to lie within local anatomical zones (Fig. 1c, right) expanded from the original ROIs to the borders of any anatomical regions that they touched. To ensure that higher connectivity of the optimized ROIs generalized to independent data, a leaveone-out cross-validation approach was adopted, with optimization performed using N Ϫ 1 animals and connectivity of the resultant ROIs measured in data from the left-out animal. Because each cross-validation fold used highly overlapping sets of animals, the final ROIs tended to be very similar across folds (Fig.  1d, right). (The initial ROI labels are not necessarily anatomically meaningful after optimization. We retain them temporarily to visualize and test the effect of optimization, but new spatial clusters will later be derived using whole-brain connectivity.) After optimization, 55 of 66 intra-MD connections increased significantly (each p Ͻ 0.05, FDR corrected) and none decreased significantly. To confirm that the optimized ROIs became better connected specifically with each other, rather than showing globally higher correlations with any region, the "networkhood" of the putative macaque MD ROIs was quantified as the mean connectivity among themselves relative to their connectivity with control regions (Fig. 1e). This measure of networkhood, tested on left-out data, increased significantly after optimization [paired t test across animals: t (34) ϭ 14.8, p ϭ 2.26 ϫ 10 Ϫ16 using six selected control regions (red bars); t (34) ϭ 14.7, p ϭ 2.64 ϫ 10 Ϫ16 using all regions outside of the expanded local zones as control regions (blue bars)]. This indicates that the refined regions indeed formed a more tightly integrated network that generalized across animals compared with the initial approximation. Whereas the IPS and sPreCS ROIs shrunk during the optimization, the other ROIs became larger. Larger regions might have stronger connectivity because of increased signal-to-noise from averaging more voxels, so one might worry that the optimization procedure could have enlarged the ROIs in a nonselective manner, rather than focusing them on the true underlying mostconnected network. We therefore assessed the networkhood of the entire expanded anatomical zones within which the MD ROIs were optimized (Fig. 1c). These expanded regions were significantly more connected then the original ROIs [paired t test across animals: t (34) ϭ 6.75, p ϭ 9.35 ϫ 10 Ϫ8 using six selected control regions (red bars); t (34) ϭ 5.96, p ϭ 9.62 ϫ 10 Ϫ7 using all regions outside the expanded local zones as control regions (blue bars)], but significantly less connected than the optimized ROIs [paired t test across animals: t (34) ϭ Ϫ6.79, p ϭ 8.30 ϫ 10 Ϫ8 using six selected control regions (red bars); t (34) ϭ 6.90, p ϭ 5.96 ϫ 10 Ϫ8 using all regions outside the expanded local zones as control regions (blue bars)]. This suggests that the optimized ROIs indeed converged on especially well connected subregions, rather than expanding indiscriminately.

Whole-brain analysis of regions connected to MD network
We next assessed the connectivity maps from the optimized macaque MD candidates to all voxels in the brain to investigate whether they showed consistent connectivity with any other regions in addition to themselves. For example, although we had focused a priori on frontoparietal ROIs, in humans these foci often coactivate with lateral occipital cortex both during visual tasks (Fedorenko et al., 2013) and at rest . Furthermore, in the macaque, it is unknown whether an MD-like network might contain more or fewer foci than are observed in humans. Again, connectivity maps for each animal were calculated using seed ROIs optimized using independent data from the other animals. Connectivity maps were then averaged across animals. Figure 2a illustrates voxelwise connectivity maps for each ROI. In addition to locally high correlations near each seed ROI, substantial similarity between the connectivity maps is apparent, as would be expected if the ROIs form part of a connected network. The overall network structure is summarized in Figure 2b, which uses a conjunction test to assess the number of MD ROIs that are significantly connected to each voxel (generalized Simes method for combining dependent p-values; Heller et al., 2007; FDR Ͻ 0.025). To derive a final definition of putative monkey MD clusters, this whole-brain conjunction map was thresholded and divided into spatially contiguous clusters in a manner similar to the partitioning of the original human MD conjunction of cognitive demand activations (see Materials and Methods). This identified eight local clusters (Fig. 2c). (Although each cluster was spatially contiguous in volumetric space, sulcal clusters could appear split across opposing banks of the sulcus when rendered onto a cortical surface.) Figure Table 1 lists the clusters and the anatomical regions with which they have the greatest overlap.
Of the eight clusters, three lie in an anterior-posterior direction along the lateral frontal surface (shades of green in Figs. 2c,  3). The most anterior of these is in the posterior principal sulcus, with a local peak in area 46. A large middle cluster covers the inferior portion of the arcuate sulcus (predominantly area 6Va and overlapping the lateral frontal eye fields). The more posterior cluster is centered on area 4C, including the arcuate spur and spreading onto the precentral and superior frontal gyri.
Another cluster (yellow in Figs. 2c, 3) is in the anterior parietal lobe, having greatest overlap with area 7b, but with some spreading onto the parietal operculum and along the intraparietal sulcus into the anterior intraparietal area and the dorsal lateral intraparietal area (LIPd). This macaque intraparietal cluster is much smaller than the intraparietal MD cluster in humans and occupies a more anterior position. The cluster partially overlaps several parietal anatomical regions, but is not especially well matched to any. It is unclear whether this might be due to methodological factors or if it reflects a genuine functional difference across species, as discussed later. The parietal cluster does encroach on region LIPd, however, where single-cell recordings show adaptive coding of task-relevant information, as expected of MD cortex (Toth and Assad, 2002;Freedman and Assad, 2006).
A further cluster (turquoise in Figs. 2c, 3) lies in the insula. This cluster covers anterior and dorsal insular cortex toward the frontal operculum, as in humans, although the macaque cluster also appears to spread more posteriorly. This corresponds well to anatomical borders, with the final cluster closely matching the anterior-posterior extent of the granular/dysgranular anatomical region.
A pair of medial frontal MD clusters (orange and red in Figs. 2c, 3), in macaques as in humans, include anterior cingulate cortex and supplementary motor areas. Although defined as a single cluster in the human, it divides into rostral and caudal clusters in the macaque, centered on areas 24d and 6M. Division into at least two macaque clusters seems appropriate given their good match to distinct anatomical regions. It should be noted that the rostral cluster covers only the anterior portion of area 24d, which some other parcellations label 24c (Matelli et al., 1991). Although the caudal cluster is continuous in volumetric space, on the cortical surface, it splits onto the banks of the sulcus, including areas 6M and 24ab, but not 24d (or 24c) lying between them in the depth of the sulcus. It is therefore unclear whether 24ab would be better grouped with the more rostral cluster or described as a separate cluster. The medial frontal clusters reported here appear to overlap the supplementary and cingulate motor areas (Paus, 2001;Dum and Strick, 2002) and the "attention-orienting" cingulate subdivision described by Hutchison et al. (2012b).
A final unexpected cluster (maroon in Figs. 2c, 3) lies along the lateral fissure, predominantly overlapping primary auditory cortex (A1), but spreading onto the parietal operculum (including the second somatosensory area, although this is noncontiguous on the cortical surface). Unlike the clusters described above, this region does not form part of the human MD system.

Comparison with previously reported anatomical projections
The CoCoMac database collates macaque axonal tract tracing data from Ͼ450 publications (Stephan et al., 2001;Bakker et al., 2012). We queried CoCoMac 2.0 to generate a list of axonal projections between all regions in the LV-FOA-PHT parcellation and to convert this to a connectivity matrix indicating connections as present, absent, or unknown. We identified the seven atlas regions that had maximum overlap with our MD functional connectivity clusters, excluding A1, and assessed the proportion of known projections involving these nodes, ignoring selfconnections and links for which connectivity was unknown. Of the MD nodes, only the insula had a known direct connection with A1. Among the MD nodes, 94% of connections were present. This compares to 55% of reported connections between these nodes and all other nodes and 51% of reported connections among all non-MD pairs. We also investigated a set of posterior parietal regions that overlapped substantially with the parietal MD region initially warped from humans, but not with the final optimized macaque clusters. A total of 57% of connections were reported between MD nodes and these posterior parietal nodes.
Comparison with anatomical connectivity thus confirmed that the putative MD clusters are physically very well connected, with the exception of A1. Although MD regions have moderate anatomical connectivity with posterior parietal regions, this is comparable to their average anatomical connectivity with all cortical regions and weaker than their connectivity with 7b and other MD nodes.

Discussion
Starting from a human frontoparietal network recruited by multiple cognitive demands, we identified corresponding, well connected regions in anesthetized macaques consistent with the proposal of a macaque counterpart to the human MD network. Because recent studies suggest evolutionary reorganization of frontoparietal cortex Patel et al., 2015), it is interesting that the broad anatomical regions of human frontoparietal MD cortex (preSMA/ACC, IPS, insula, and multiple lateral prefrontal foci) are also all consistently connected in anesthetized macaques. In addition, apparent interspecies differences in the relative size and location of these connected clusters are noteworthy.
The proposed macaque MD system is compatible with suggestions of evolutionarily novel human frontoparietal networks. Mantini et al. (2013) partitioned human and macaque restingstate networks using independent component analysis and then performed cross-species matching using spatial correlation, leav-   ing three human frontoparietal and cingulate components unmatched. The putative macaque MD network appears to include several human-macaque homologs identified by  to match human "dorsal attention," "lateral prefrontal," and "language" components. The present cingulate, posterior lateral frontal, and anterior intraparietal clusters are not well covered by the macaque resting-state networks selected by . Other studies suggest a generally well conserved organization of frontal and parietal cortex (Mars et al., 2011;. In humans, the MD system may be supplemented by either novel frontoparietal networks or stronger connections that might reflect evolutionary expansion of prefrontal and parietal cortex (Van Essen and Dierker, 2007). The macaque intraparietal cluster was surprisingly small and anterior compared with its human counterpart. In anesthetized macaques, it is not uncommon to observe relatively weak connectivity of posterior intraparietal cortex with many frontal regions (Hutchison and Everling, 2013) or weaker connectivity than expected from human data (Mars et al., 2011). Potentially, posterior parietal cortex may be particularly susceptible to reduced connectivity under anesthesia (Schröter et al., 2012) and resting-state networks spanning prefrontal and posterior intraparietal cortex may be especially apparent in awake animals . Alternatively, the absence of posterior intraparietal cortex among proposed MD regions may reflect meaningful interspecies differences in posterior intraparietal function (Orban et al., 2004). Supporting this possibility, sus-tained attention activates multiple human parietal foci, but a single anterior focus in macaques (Patel et al., 2015). Similarly, interspecies activity correlations of macaque posterior parietal regions (e.g., parietal occipital, medial intraparietal, and posterior intraparietal areas) match human occipital regions better than human parietal regions .
Structural connectivity likely constrains covariance between time series, especially at long time scales (Honey et al., 2007), and moderate correlations between structural and functional connectivity have been reported (Miranda-Dominguez et al., 2014). We therefore considered how the present results correspond to anatomical connections, confirming that the putative frontoparietal MD clusters are well connected physically as well as functionally. They are not anatomically well connected to A1, suggesting that the unexpected cluster in the lateral fissure may not belong to the core network. They have only moderate direct projections with more posterior parietal regions, suggesting that the unexpected absence of these regions among the functional connectivity clusters has some anatomical support and may not be merely an effect of anesthesia. Indeed, although human and macaque structural connectomes are broadly similar, parietal regions have relatively large differences in their connectivity patterns (Goulas et al., 2014).
Cortical subdivisions, whether structurally or functionally derived, can be described hierarchically, with regions grouped or split depending on hierarchical level (Modha and Singh, 2010; Hutchison and Everling, 2013). The MD network represents an intermediate level of description that is more specific than the densely connected core that spans task-positive and task-negative regions (Modha and Singh, 2010;Ercsey-Ravasz et al., 2013), yet encompassing many cytoarchitectonically distinct areas (Van Essen et al., 2012) and with finer-grained functional structure evident in humans (Dosenbach et al., 2007) and macaques (Wardak et al., 2010;Mars et al., 2011;Hutchison et al., 2012b;Sallet et al., 2013;. Given the obvious spatial peaks within human MD cortex, we have suggested a similar spatial subdivision of macaque MD cortex. Alternative subdivisions are possible, however, and different levels of description are likely informative in different contexts, such as predicting function or disease. Whole-brain connectivity of MD ROIs revealed one nonfrontoparietal cluster in the lateral fissure overlapping primary auditory cortex. In humans, this region is neither activated by diverse cognitive demands (Fedorenko et al., 2013) nor well connected to executive control networks . Conversely, human ventrolateral occipitotemporal cortex commonly responds to cognitive demand (Fedorenko et al., 2013) and correlates with executive control networks at rest , yet was neither strongly nor consistently connected with macaque MD ROIs. One might speculate that visual stimuli dominate human experiments and visual imagery may occupy human thought at rest, producing coupling between human MD regions and occipitotemporal visual representations. For lightly anesthetized macaques, scanner noise may instead be the dominant stimulus being processed as awareness fluctuates, driving coupling between macaque MD and auditory cortex. One human study reported that sedation partially preserved frontoparietal connec- Table 1. Description of regions with significant connectivity to five or more optimized MD ROIs, clustered by local peaks of mean connectivity as shown in Figure 2c and Figure 3 Location of cluster  Region 24ab is not contiguous with the rest of the cluster on the cortical surface, so could alternatively be grouped with the rostral anterior cingulate cluster or described as a separate cluster.
tivity of executive control networks, whereas occipitotemporal connectivity from a right hemisphere executive control seed was no longer detected (Boveroux et al., 2010). In macaques, auditory cortex has been found to correlate with lateral and dorsal prefrontal cortex under anesthesia (Hutchison and Everling, 2013), whereas little auditory-with-lateral-prefrontal functional connectivity is observed in awake animals . Given the numerous challenges of fMRI in behaving macaques, few studies specifically contrast levels of cognitive demand. In humans, inhibiting a prepotent response produces typical MD activation. In macaques, antisaccades (away from targets that normally attract gaze) evoke fMRI responses consistent with the proposed MD regions, including the principal, arcuate, anterior cingulate, and intraparietal sulci (Ford et al., 2009). An earlier study investigating response inhibition during a Wisconsin Card Sorting Test also observed activations partially overlapping the current MD clusters, including the inferior arcuate sulcus and insula (Nakahara et al., 2002). Whereas those studies compared conditions differing objectively in difficulty, others observe frontoparietal activation while carefully matching behavioral performance. These implicate frontoparietal regions resembling the current MD clusters in aspects of attentional control: searching for salient targets (Wardak et al., 2010) and covertly shifting attention (Caspari et al., 2015). This matches human frontoparietal involvement in attention switching (Greenberg et al., 2010) and detecting changes to attended items (Beck et al., 2001;Hon et al., 2006;Finoia et al., 2015), which often covary with cognitive demand, but also produce frontoparietal activation when performance is equated (Greenberg et al., 2010) or irrelevant (Hon et al., 2006; for a macaque equivalent, see Stoewer et al., 2010).
Given the restrictions on human tracer studies and the challenges of fMRI in behaving animals, functional connectivity provides an important bridge between human and animal models (Hutchison and Everling, 2012;Miranda-Dominguez et al., 2014). Monkeys, unlike humans, are usually anesthetized, restrained, extensively trained, or rewarded. Variability in cognitive state therefore confounds interspecies comparisons and requires better understanding. Despite consensus that functional networks are largely preserved under anesthesia in humans, macaques, and rodents, albeit with weakened connectivity Greicius et al., 2008;Boveroux et al., 2010;Mars et al., 2011;Hutchison et al., 2012a;Liang et al., 2012;Sallet et al., 2013;Miranda-Dominguez et al., 2014), the extent to which arousal has differential effects across the brain or upon short-range versus long-range connections is less certain (Moeller et al., 2009;Boveroux et al., 2010;Liang et al., 2012;Schröter et al., 2012). Similarly, although analyses attempt to control rigorously for physiological confounds, restraints, sedatives, contrast agents, and external ventilation might have residual effects. Finally, differences in the reliability and spatial distribution of the fMRI signal (e.g., from magnetic field inhomogeneities or surface coil positioning) may affect interspecies comparisons.
More fundamentally, human MD regions are defined from contrasts of cognitive demand and, ultimately, we want to understand both human and macaque brain systems in terms of their role in cognition. Comparing MD clusters with available macaque fMRI studies of response inhibition and attention switching informally suggests similar functional properties to their human counterparts. Full validation of the proposed macaque MD network, however, requires fMRI studies that test explicitly these regions' specificity in responding to diverse cognitive de-mands. The putative macaque MD network is made available for this purpose (http://imaging.mrc-cbu.cam.ac.uk/imaging/ MDsystem). Having delineated the MD network, simultaneous electrophysiological recordings from its multiple subregions can investigate coordination of activity at finer temporal scales in the service of goal-directed behavior.
Overall, functional connectivity suggests that a macaque homolog of the human MD system exists, even under anesthesia and despite interspecies differences in frontoparietal expansion (Van Essen and Dierker, 2007), task responses (Patel et al., 2015), and some resting-state networks . Beyond advancing primate comparative neuroscience, the putative macaque MD regions can guide electrophysiological recordings or be used in task-based fMRI to test predictions of similar functional properties to human MD cortex. We hope that characterizing the relationship between human and macaque MD systems will enhance the translational potential of nonhuman primate models in understanding human brain function and pathology.

Notes
Supplemental material for this article is available at http://imaging.mrccbu.cam.ac.uk/imaging/MDsystem. The putative macacque frontoparietal MD regions, their conjunction connectivity map, and the final connected clusters, are available to download. ROIs are provided after initially warping from humans and after subsequent optimization. Data can be downloaded in the space of the 112RM-SL template and the F99 brain volumes. Files are also provided for visualizing the data on the F99 cortical surface using Caret software (http://brainvis.wustl.edu/). This material has not been peer reviewed.