Structural mediation of human brain activity revealed by white-matter interpolation of fMRI
Anjali Tarun, Hamid Behjat, David Abramian, Dimitri Van De Ville
SStructural mediation of human brain activity revealedby white-matter interpolation of fMRI
Anjali Tarun a,b,1 , Hamid Behjat c , David Abramian d,e , Dimitri Van DeVille a,b a Institute of Bioengineering, ´Ecole Polytechnique F´ed´erale de Lausanne, Geneva, 1202,Switzerland b Department of Radiology and Medical Informatics, University of Geneva, Geneva, 1202,Switzerland c Center for Medical Image Science and Visualization, University of Link¨oping,Link¨oping, 58183, Sweden d Department of Biomedical Engineering, University of Link¨oping, Link¨oping, 58183,Sweden e Department of Biomedical Engineering, Lund University, Lund, 22100, Sweden
Abstract
Anatomy of the human brain constrains the formation of large-scale func-tional networks. Here, given measured brain activity in gray matter, weinterpolate these functional signals into the white matter on a structurally-informed high-resolution voxel-level brain grid. The interpolated volumesreflect the underlying anatomical information, revealing white matter struc-tures that mediate functional signal flow between temporally coherent graymatter regions. Functional connectivity analyses of the interpolated vol-umes reveal an enriched picture of the default mode network (DMN) andits subcomponents, including how white matter bundles support their for-mation, thus transcending currently known spatial patterns that are limitedwithin the gray matter only. These subcomponents have distinct structure-function patterns, each of which are differentially recruited during tasks,demonstrating plausible structural mechanisms for functional switching be-tween task-positive and -negative components. This work opens new avenuesfor integration of brain structure and function, and demonstrates how globalpatterns of activity arise from a collective interplay of signal propagationalong different white matter pathways.
Email address: [email protected] (Anjali Tarun) a r X i v : . [ q - b i o . N C ] A ug ntroduction Coordination of distant neuronal populations gives rise to a vast reper-toire of functional networks that underpin human brain function. Usingfunctional magnetic resonance imaging (fMRI), temporally coherent activitycan be investigated using measures of functional connectivity (FC). On theother hand, the mediation of inter-regional communication by the anatomi-cal scaffold can be conveniently summarized by structural connectivity (SC)extracted from diffusion-weighted MRI (DW-MRI). Over the past decade,several methods have been proposed to bridge the gap between SC and FC.Seminal works have found evidence for strong statistical interdependence be-tween separately defined SC and FC [1–4]. Limited by the bivariate andsummarizing nature of the analysis, the effect is capturing only a generaltrend of correlation. Following after were studies that specify regions of in-terests (ROI) as a FC prior for extracting white matter pathways [5, 6], mostof which were specifically concentrated to the analysis of the default modenetwork (DMN), a set of brain regions that are known be more engaged dur-ing rest [7]. In contrast to extracting SC from FC priors, a number of studieshave attempted to reproduce brain activity from predefined structural con-nectomes through numerical simulations [8–13].Studies that extract SC from FC or vice versa are mostly hypothesis-driven and entail many explicit assumptions. In order to understand how dis-tributed patterns of functional activity arise from a fixed underlying anatomy,a need for research methodologies that are observer-independent and data-driven are of utmost importance. A common approach for data-driven ap-proaches are based on blind-decomposition techniques. By combining dif-fusion anisotropy ( e.g. , fractional anisotropy, axial and radial diffusivities)and classical FC to build a joint SC-FC measure [14, 15], structural andfunctional alterations in healthy and clinical populations are extracted us-ing multimodal canonical correlation analysis and joint independent compo-nent analysis (ICA). On the same line, using ICA, the extraction of concur-rent white matter (WM) bundles and gray matter (GM) networks from atractography-based data matrix that records the number of streamlines run-ning from a GM position to a WM location has been proposed [16]. Alsousing tractography, Calamante and colleagues [17] introduced the concept oftrack-weighted dynamic FC (TW-dFC), where as a first step, WM stream-2ines from two structurally connected voxel end-points are obtained. Then,these streamlines are weighted according to the FC at their endpoints over asliding-window that is arbitrarily assigned. TW-dFC then produces a set offour-dimensional volumes showing the averaged FC between correspondingstreamline end points onto the WM at the temporal resolution of the win-dow used. While these studies have made considerable progress in linkingstructure and function, none of them jointly and simultaneously integratesfunctional and diffusion data into a single integrated framework, which al-lows the natural emergence of full-brain spatial patterns that covers both theWM and the GM.Lately, to further transcend our current understanding of structure-functionrelationships, attention has been set on studying SC and FC through the lensof more technical frameworks borrowed from other research fields, such aspropagator-based methods [18–20], and control network theoretical tools [21].Graph signal processing (GSP) for neuroimaging is another emerging field [22],where initially, a graph is defined by identifying regions in GM as nodes, andmapping strength of their SC through WM to edge weights. Functional dataare then interpreted as time-dependent graph signals, on which connectome-informed signal processing operations can be performed. This frameworkincorporates connectivity through WM, but only as SC between pairs of GMregions.In this work, we advance the GSP concept by modeling WM explicitlyas nodes of the graph for which local connectivity is known from DW-MRI.This allows relating measures of brain activity on the GM with their media-tion through WM; i.e., how particular patterns of brain activity jointly relyon an ensemble of WM pathways. This is done by extending the existingapproach from region-wise analysis to a whole-brain, voxel-level perspective.We define the blood–oxygenation level dependent (BOLD) time-series fromresting-state fMRI as dynamic signals residing on the nodes of a high resolu-tion (i.e., 850,000 voxels) brain graph. By generalizing principles of classicalsignal processing in regular domains [23–26] to irregular graphs [27, 28], weinterpolate functional signals, measured on the GM, into the WM, usinga whole-brain voxel-wise connectome to guide the process. The recoveredvolumes provide a novel perspective on how coherent GM regions commu-nicate through various combinations of WM pathways. The application ofconventional and dynamic FC tools on the recovered volumes reveals newstructure-function relationships of DMN-related networks, illustrating thestructural mechanism for the dynamic switching between task-positive and3negative subsystems of the DMN in different phases of working memory andrelational task paradigms.
Results
Recovery of activity in white matter that is compatible with structural con-nectivity
Signal recovery is a classical task that relates to image denoising [23], sig-nal inpainting [24, 25], and sensing [26]. Here we embed the three-dimensionalvoxel-level grid into a graph where nodes are voxels, and edges of the 26-neighborhood of each node encode strength of structural connectivity asmeasured by DW-MRI. Functional volumes are preprocessed and upsampledto the same resolution as this structural grid. Nodes in the GM are assignedthe values from these functional volumes, while the remaining nodes, in par-ticular, those in the WM, remain unassigned. The values of the unassignednodes are then recovered by solving an optimization problem that relies ontwo assumptions: (i) values in the GM nodes should strongly influence thedata fitness (and thus remain reasonably unchanged); (ii) whole-brain sig-nals should maintain smoothness according to the graph structure. Theassumption of smoothness with respect to the graph ensures that signals areinterpolated according to the brain’s structural scaffold. The overall pipelineis illustrated in Fig. 1(A) and explained in details in the Supplementary In-formation. For an illustrative frame of resting-state fMRI, we observe strongactivity in the posterior cingulate cortex (PCC) and medial prefrontal cortex(mPFC) for the original and interpolated volumes (Fig. 1 (B)). The latter re-veals additional WM signals showing high activation in the cingulum bundle(R/Lcing), atlased using NeuroImaging Tools and Resources Collaboratory(NITRC) repository. In addition to major WM bundles that are clearly cap-tured, we also find short-range WM connectivity within the cortical foldingsof the GM [29]. The original PCC-averaged signal (Fig. 1(C)) shows closesimilarity to its interpolated counterpart, a direct consequence of the secondconstraint imposed in the signal recovery formulation ( i.e. , retaining the sig-nals within the GM nodes fixed). In contrast, we find a huge difference in themean signals within the R/Lcing. PCC, mPFC and recovered R/Lcing sig-nals all show similar trends, while the original R/Lcing signal behaves muchmore randomly. This observation demonstrates that interpolated functionalsignals are smooth over the structural brain grid, reflecting in this particu-4 riginal BOLD Masked GM BOLD interpolated BOLDDiffusion Spectrum Imaging (DSI)GM BOLD imposedon whole brain grid interpolated PCCraw PCC r a w / i n t e r po l a t ed B O L D time raw mPFCinterpolated cingulumraw cingulum BC BOLD signal A interpolated BOLD signal PCC cingulummPFC +3.0 -3.0 S . D . z = 23z = 23 y = -28y = -28 x = -10x = -10 Figure 1:
Workflow of the graph signal recovery framework. (A) GM BOLDsignals are extracted from fMRI volumes, one signal per time instance, through maskingthe volume by the thresholded probabilistic GM map (threshold=0.3). ODFs associated toall voxels are extracted from the diffusion MRI data. The ODFs are then embedded on athree-dimensional, 26-neighborhood connectivity grid, forming a probabilistic connectomeat voxel-level resolution. A signal value can be associated to each voxel on the grid. Toinitiate signal recovery at each time instance, the associated GM BOLD signal is mapped tothe GM voxels within the grid, and the remaining voxels are set to zero. The interpolatedBOLD volume is obtained through minimizing a cost function that balances between (1)retaining the GM signal fixed and (2) obtaining a smooth signal over the entire brain gridrelative to the underlying structure. (B) Axial, sagittal, and coronal views of a singleslice of a BOLD volume is compared to its corresponding interpolated BOLD volume.(C) Original and interpolated BOLD signal time-series of PCC, mPFC and cingulum; theselected time frame in (B) corresponds to the highlighted time point.
PCC-seed connectivity map
One conventional method to explore FC is seed-correlation analysis, wherethe BOLD time course of an ROI is correlated to the activity profiles of allbrain voxels to locate temporally coherent areas. In order to capture thejoint structural-functional connection of the DMN, we perform a PCC seed-connectivity analysis [30] on the interpolated volumes derived from four ses-sions of resting-state scans in a population of 51 subjects (244,800 volumesin total). The superior axial slice ( z =25) in Fig. 2(A) and (B) shows thePCC and the mPFC, including the bilateral inferior parietal cortex (IPC).In addition to the expected GM FC pattern, Fig. 2(B) also shows the distinctWM structures that support the structural connections between these tem-porally coherent cortical regions. The most prominent connections are theR/Lcing that connect the PCC and the mPFC, the forceps minor (Fminor)that provides intra-connectivity between gyri within the mPFC, and the leftand right superior longitudinal fasciculi (R/Lslf) that support the long-rangeconnection between the IPC, and the posterior and frontal regions. Interest-ingly, we also find the two-sided corticospinal tracts (R/Lcst) that traversethe brainstem all the way to the motor cortex, the genu of the corpus cal-losum connecting the two hemispheres, as well as the bilateral hippocampalcingulum (R/Lcing2) that exit the caudal aspect of the PCC and continuetowards the medial temporal lobe (MTL). Structural mediation of dynamic functional connectivity
The volume shown in Fig. 1(B) corresponds to a time point with signifi-cant activity in the PCC (green shade, Fig. 1(C)). FC obtained from conven-tional seed-correlation analysis can be approximated by averaging all framesfor which the seed has high activity [31]. This congregation of significantlyactive frames can also be temporally clustered into distinct and meaning-ful co-activation patterns (CAPs). Using this approach, we decomposed theDMN pattern into eight PCC-related co-activation patterns (CAPs). Fig. 3presents the CAPs obtained from the interpolated volumes, termed inter-polated CAPs (in-CAPs), as well as the CAPs generated from the originalBOLD signal, denoted GM CAPs. The resulting spatial patterns reveal a setof WM structures unique to each in-CAP. In-CAP 1 is the most frequently6 = - 32 z = 10 z = 25 y = - 9y = - 21y = -47x = - 21 y = -33 AB BOLDinterpolatedBOLD y = 0 pe r c en t il e2 pe r c en t il e Lcing2 LcingLcst LslfLilfLslf corpus callosum corpus callosumFminor x = - 32 z = 10 z = 25 y = - 9x = - 21 y = -33
IPC mPFC PCC C MTL
Figure 2:
PCC-seed correlation for the (A) original BOLD volumes and the (B) interpo-lated volumes, averaged across all subjects. WM bundles connecting temporally-coherentGM regions are uncovered from interpolated volumes. (C) 3D view of the observed WMstructure.
234 5678 +1.5-1.5
S.D -21 -9 3 15 27 39 51 -21 -9 3 15 27 39 5123.8%
Occurrences (%)CAP Number G M C A P i n - C A P G M C A P i n - C A P G M C A P i n - C A P G M C A P i n - C A P G M C A P i n - C A P G M C A P i n - C A P G M C A P i n - C A P G M C A P i n - C A P Figure 3:
PCC-seed co-activation patterns (CAPs).
Eight co-activation maps ob-tained by temporal decomposition of the top 15% significant frames related to the posteriorcingulate cortex (PCC), extracted from (1) the original BOLD volumes (GM CAPs) and(2) interpolated volumes (in-CAPs). The in-CAPs are numbered according to their fre-quency of occurrence. Visual inspection of the in-CAPs reveals strong GM similarity withthe observed GM CAPs. We found eight distinguishable structure-function networks re-lated to the default mode network (DMN), each of which varies in terms of the observedWM structures that conjoin distinct PCC-coherent GM regions.
Relevance of WM interpolated patterns during rest and task
Whereas GM CAPs reveal instantaneous co-activation of multiple brainregions, in-CAPs complete the structure-function picture through the ad-dition of distinct and functionally relevant WM structures that interposespatially-distant GM co-activations. The concurrence of the two CAP typesis confirmed by their generally matching frame assignments (Supplemen-tary Fig. S1). This strong temporal correspondence illustrates the abilityof the proposed method in successfully capturing distinct WM structuresconnecting temporally varying PCC-related networks, without compromis-ing information from the original data. To further study in-CAPs and theircorrespondence with the DMN subsystems, we computed the average signalwithin major WM tracks (Fig. 4(A)) using the NITRC WM atlas, and com-pared them to the main components of the DMN (Fig. 4(B)), as classifiedby Andrews-Hannah and colleagues [32]. We observed good correspondencein the spatial patterns of GM regions and WM tracks, especially in terms ofactivation polarity from the dorsal (positive-valued) to the ventral (negative-valued) side. In-CAP 8, however, is positive across all DMN sub-components,and in contrast to all other in-CAPs, also shows a positive bilateral uncinatefasciculus (L/Runc).Previous work have found that activations in both the hippocampusand ventral mPFC during retrieval-mediated learning support novel infer-ence [33]. We hypothesize that the positive activity of L/Runc in in-CAP9 reflects its potential role in mediating the coherent activity between thepositive temporal poles, hippocampus, and ventral mPFC. To verify thisassociation and to further explore the functional meaning of the obtained in-CAPs, we computed their occurrence probability upon working memory andrelational tasks (Fig. 4(C)). To do so, we interpolated signals into the WMas before, extracted frames with high activation in the PCC, and assignedthem to their closest in-CAP. We found that although in-CAP 1 is the mostrecruited at rest, it mostly occurs during task blocks (as opposed to rest in-tervals). In contrast, in-CAP 8 occurs the least during resting state, but it isvividly recruited during the rest epochs of the task paradigms. These resultsare evaluated using a two-factor ANOVA on the occurrences of the differentinterpolated CAPs during rest and tasks. Results showed highly significantmain effects of CAP type and task, as well as the interaction between both,with all corresponding p-values less than 0 . Discussion
General findings
We modeled the brain grid, together with local, voxel-to-voxel probabilis-tic connections based on DW-MRI, as a large graph onto which brain activityin GM is interpolated into WM. The interpolated brain volumes remarkablyreflected WM structures that support distributed patterns of activity in GM.We explored the neuroscientific relevance of the interpolated volumes by ap-plying conventional and dynamic FC tools, namely seed-correlation [30] andCAP analysis [31], using the PCC as a seed. The resulting in-CAPs unrav-eled the exquisite complexity of recruited WM patterns underlying typicallyobserved GM co-activations (Fig. 3). In-CAPs also showed strong neuro-physiological relevance, as they revealed the segregation of the DMN intotask-positive and task-negative sub-components (Fig. 4C).
Structure-function relationships of PCC-related networks
Interpolating functional signals measured on GM into the WM has re-vealed new structural-functional relationships that are akin to the DMN.Fig. 4(A) shows the signal average within WM bundles observed in ourPCC seed-correlation map, consistent with the visually observed tracks inFig. 2. Prominent tracks include the R/Lcing, the R/Lslf, and the Fminorin close agreement with previous findings of WM structures associated with10he DMN [5, 6]. Although seed-correlation analysis is a straightforward mea-sure of FC, it is debased by its transitive nature. The connectivity obtainedfrom this static approach, therefore, summarizes all major WM connectionslinked to the PCC. On the other hand, CAP analysis can account for time-dependent behavior; thus, the in-CAPs reveal different sets of WM bundlesthat are characteristic of each of the observed GM CAPs at different time-points.While it has been established that the DMN anatomically consists of theanterior and posterior midline, the bilateral parietal cortex, prefrontal cortex,and the MTL [34], these regions have been found to functionally dissociateaccording to the ongoing internal processes [32]. We surmise that this func-tional disintegration is captured by in-CAPs, given that distinct in-CAPsoccur in different phases of working and relational memory tasks (Fig. 4(C)).The MTL activates when internal decisions involve constructing a mentalscene based on memory [32]. However, despite being dominantly active inthis region, in-CAP 8 appears to be more active during blocks of rest, and in-stead, in-CAPs 1 and 2 dominate in periods of task blocks when the stimulusis presented. These findings suggest that in-CAP 8 is a characteristic task-negative structure-function network that activates during the maintenancephase ( i.e. , when subjects are focused on consolidating an observed stim-ulus), and conversely, in-CAPs 1 and 2 are task-positive structure-functionnetworks that are engaged during the encoding and retrieval periods whenthe external stimulus is presented. Altogether, our results provide new in-sights into the well-established dichotomy of the human brain [35, 36], andclearly demonstrate that the DMN and task-positive networks ( e.g. , attentionnetwork, working memory network) are not exclusively expressing oppositeactivity [37, 38].Arranging the WM tracks and the DMN sub-components according tothe polarity of their activity in each of the in-CAPs (Fig. 4(A-B)) reveals ageneral trend of positive dorsal areas, such as those observed from Fmajor toLilf, and from left PCC (l-PCC) to left temporal parietal junction (l-TPJ).On the other hand, ventral components (e.g., Rcing2 to Runc, left lateraltemporal cortex (l-LTC) to left temporal pole (l-Temp)) are more varied,showing a subtle divide between more negative ventral (in-CAPs 1-4) versusmore positive ventral default networks (in-CAPs 5-8). Sagittal slices (Sup-plementary Fig. S3) show strong activation in the ventral PCC of in-CAP 8,together with the R/Lcing2, which runs from the caudal aspect of the retro-splenial cortex (RSP) to the hippocampal formation (HF). The bilateral HF,11 major FminorFornix LcingLcing 2
Lcst
Lifo LilfLslf Lunc A in-CAP occurrence in rest and tasksSpatial Average within WM bundles F m a j o r F m i no r L c i ng L s l f R s l f R cs t L cs t R c i ng R il f L il f F o r n i x R c i ng2 R i f o L c i ng2 L i f o Lun c R un c seed corrin-CAP 1in-CAP 2in-CAP 3in-CAP 4in-CAP 5in-CAP 6in-CAP 7in-CAP 8 -1.5-1.0-0.500.51.01.5 B Spatial Average within GM DMN sub-components l - P C C r- P C C r- T P J l - R s p l - p I P L r- R s pd M P F C l - a M P F C r- a M P F C r- p I P L l - T P J l - L T C r- L T C r- P H C l - P H C r- T e m p P r- H F l - H F v M P F C l - T e m p P in-CAP 1in-CAP 2in-CAP 3in-CAP 4in-CAP 5in-CAP 6in-CAP 7in-CAP 8 -2.5-2-1.5-1-0.500.511.522.5 C in-CAPs P r babab ili t y d i s t r i bu t i on PHC HF vMPFC aMPFCdMPFCpIPLTempP
RestWorking MemoryRelational £21,200 TR TR N/AR
Figure 4:
Spatial and temporal characteristics of in-CAPs.
The overall spatialpattern of the in-CAPs is explored through spatial averaging of signals within (A) 17 majorWM bundles, and (B) GM sub-components of the DMN [32]. The arrangement of WMstructures and DMN sub-components based on their polarity show a concordant patternof positive-to-negative transition from dorsal-to-ventral regions. (C) The probability ofoccurrence of in-CAPs in resting-state, working memory, and relational tasks. Each barin the plot corresponds to the probability of that particular in-CAP to appear in the top15% significant volumes of rest and task scans. The stack of colors denotes the proportionof frames that occur during blocks of rest (R) when there is no stimulus presented, andduring blocks of task (T) in which subjects are expected to perform memory or relationaltest.
Methodological perspectives
The proposed methodology is one of the very few studies that integratesand jointly models functional and diffusion data in an observer-independentand fully data-driven manner. Previous works have relied on blind-decompositiontechniques and the use of either diffusion anisotropy measures or tractogra-phy algorithms to approximate SC [14–16]. However, ICA’s principal as-sumption of spatial independence leads to components that have low spa-tial similarity, while measures of diffusion anisotropy are unable to resolvecrossing fibers as they are extracted from low-order diffusion tensor mod-els, making their interpretation misleading and prone to erroneous conclu-sions [42]. Additionally, the high susceptibility of tractography algorithmsin generating false positive [43] suggests cautious to be taken in interpretingtractography-based results. In contrast, our proposed method does not relyon tractography thereby, not only bypassing the complicated task of param-eter optimization for extracting WM tracks but also not requiring a parcel-lation scheme to define ROIs. Our approach keeps the analysis close to theDW-MRI data by working directly with local reconstructions of orientationdistribution functions, which are consequently encoded, both in direction andmagnitude, to build the brain graph. Furthermore, it is worth noting thatdespite the superficial similarity, there are fundamental differences betweenour approach and that of the TW-dFC [17]. Specifically, both approachesattempt to give a measure of FC into the WM using GM as a constraint.However, TW-dFC relies on two successive use of classical methods (e.g.,tractoragphy and sliding-window FC), while our work introduces an entirelynew concept where recovery of mediation by SC is incorporating the FCitself in a single integrated framework; i.e., the complete pattern of brainactivity in GM is constraining the recovery of WM patterns. In addition,our method works at the single-frame level and does not require to choose atemporal window to render the approach dynamic.13 ew research avenues for structure-function studies
GSP approaches as applied to functional neuroimaging rely on a graphshift operator (typically the Laplacian or the adjacency matrix) that encodesthe brain’s anatomical information. Its eigendecomposition produces an or-thogonal set of eigenvectors, termed eigenmodes , some of which are reminis-cent of well-established functional networks [20] and are useful in effectivelyestimating the strength of inter-hemispheric interactions in the brain [19].Moreover, efficient anatomically-informed decompositions of fMRI data us-ing tractography-based structural connectome [20], as well as topology en-coding GM graphs [44] have been proposed. Brain eigenmodes can be viewedas basic building blocks of increasing spatial variation along the structuralbrain graph, akin to sinusoids in classical Fourier analysis. Here, we extendthe traditionally region-level eigenmodes (defined on a limited subset of GMnodes) to a whole-brain (GM and WM), voxel resolution setting (see Supple-mentary Fig. S4), thus enabling to reconstruct structure-function networksat an unparalleled spatial resolution. As an increasing number of operationsare generalized from classical signal processing to the graph setting [22, 27],promising avenues arise to explore the many facets of brain structure andfunction.We therefore foresee two direct avenues for future research. The firstavenue aims on leveraging the proposed interpolated fMRI data. The inter-polated volumes entail additional informative voxels that extend beyond theGM, which is in particular interesting in light of the limitations and skepti-cisms in interpreting WM BOLD data [45, 46]. Previous works have demon-strated the possibility to capture functionally relevant information from theWM BOLD [47, 48], despite the well-established findings on the differencesof hemodynamic responses in GM and WM [49, 50]. At its current form,we interpolate functional signals into the WM as solely constrained by thefunctional signals from the GM. Alternatively, one may rather consider com-bining signal interpolation with signal recovery of weak signals in the WM.In our current application, the goal is to observe WM pathways that supportdistributed patterns of FC in GM, and not on the functional organizationof WM itself. The interpolated volumes can then be readily explored usingexisting tools for dynamic FC analyses, such as sliding-window correlation,ICA [51] and principal component analysis [52], to probe functional braindynamics at a much larger scale. The second foreseen research direction isto exploit the proposed anatomically-informed brain grid to implement novelGSP operations on fMRI data. High-resolution eigenmodes enable spectral14raph-based analysis of task-based and resting-state fMRI data at an un-precedented level of detail, and across the whole brain, in contrast to theconventional region-wise analysis, which is typically also limited to the corti-cal layer. Anatomically-informed spectral graph decomposition of fMRI datausing these high-resolution brain eigenmodes is anticipated to open up newperspectives on the brain’s functional organization not only within GM butalso WM. As a case in point, FC has often been associated to Euclideandistance [53, 54]; thanks to our framework, the notion of distance becomesmore meaningfully defined in terms of the spectral representation of func-tional signals [55], enabling a deepened understanding of the roles of short-and long-range connectivity in improving the efficiency of inter-areal braincommunication.In conclusion, this study presents a new framework for studying structure-function relationship that has several key advantages. The interpolationof fMRI activity into the white matter enables observing key WM struc-tures that link interacting GM regions. We believe that the introduction ofhighly-resolved human brain eigenmodes can (i) shift the existing trend ofconstructing region-wise connectomes to that of voxel-vise connectomes and(ii) expand the use of GSP repertoire in the context of functional brain imag-ing. More importantly, we anticipate the proposed joint structure-functioncharacterization to offer unprecedented benefits for the study of clinical pop-ulations hindered by complex mixes of structural and functional alterations,such as in patients with agenesis of the corpus callosum or virtual calloso-tomies [56].
Methods
Data
We used a total of 51 subjects obtained from the publicly available HumanConnectome Project dataset (HCP 1200 release), WU-Minn Consortium.Preprocessed diffusion MRI, original functional MRI, and anatomical dataare downloaded following a random subject selection scheme while maintain-ing a uniform demographic distribution (male and female, ages 22-50). MRIacquisition protocols of the HCP and preprocessing guidelines for diffusionMRI are fully described and discussed elsewhere [57]. No further preprocess-ing was applied to the diffusion data. We then extracted the ODFs usingDSI studio (http://dsi-studio.labsolver.org). On the other hand, original and15nprocessed functional MRI sequences underwent a standard preprocessingprotocol [58].
Building the voxel-wise brain grid
Similar to classical brain graphs, we define our grid as G := ( V , A ),where V = { , , , ..., N } is the set of N nodes representing the brain voxels(i.e., N = 700 −
900 thousand), and A ∈ N × N is an adjacency matrixthat encodes the likelihood of water molecules to diffuse from their currentposition to neighboring voxels. We used an ODF-based weighting-schemeon a three-dimensional 26-neighborhood mesh. The topology of the braingrid is validated by defining the graph Laplacian matrix in its symmetricnormalized form ( L = I − D − AD − ), whose eigendecomposition leads toa complete set of orthonormal eigenvectors that span the graph spectraldomain with their corresponding real, non-negative eigenvalues. The top 11eigenfunctions corresponding to the lowest spatial variation are provided inSupplementary Fig. S4. Graph signal recovery formulation
We consider a cost function that includes a least squares data-fitting termequal to the residual sum of squares (RSS), and an L2 regularization termthat takes into account smoothness with respect to the brain grid, given by˜ x = arg min x (cid:107) y − Bx (cid:107) + λ x T Lx , (1)where y is a vector of length N containing initial BOLD values within GMnodes and zero elsewhere, B is an indicator matrix that selects the GMnodes, and λ is the regularization parameter. L is the graph Laplacian,which encodes the voxel-level anatomical structure. The cost function bal-ances between (1) minimizing the RSS and retaining the original GM signals,and (2) imposing smoothness with respect to the structure of the graph. Thebalance is dictated by the optimal λ (see Supplementary Fig. S5 for detailson deriving an optimum value). The solution to Equation 1 gives the inter-polated volume for one time-point. The interpolation is done repeatedly foreach volume in the functional data. References [1] Honey, C. J. et al.
Predicting human resting-state functional connectiv-ity from structural connectivity.
Proceedings of the National Academy of ciences of the United States of America , 2035–2040 (2009). URL .[2] Andrews-Hanna, J. R. et al. Disruption of Large-Scale Brain Systemsin Advanced Aging.
Neuron , 924–935 (2007).[3] Horn, A., Ostwald, D., Reisert, M. & Blankenburg, F. The structural-functional connectome and the default mode network of the humanbrain. NeuroImage , 142–151 (2014).[4] Supekar, K. et al.
Development of functional and structural connectivitywithin the default mode network in young children.
NeuroImage ,290–301 (2010).[5] Greicius, M. D., Supekar, K., Menon, V. & Dougherty, R. F. Resting-state functional connectivity reflects structural connectivity in the de-fault mode network. Cerebral cortex , 72–78 (2009).[6] Van Den Heuvel, M., Mandl, C. W., Kahn, R. S. & Hulshoff Pol, H. E.Functionally linked resting-state networks reflect the underlying struc-tural connectivity architecture of the human brain. Human Brain Map-ping , 3127–3141 (2009).[7] Greicius, M., Krasnow, B., Reiss, A. & Menon, V. Functional con-nectivity in the resting brain: a network analysis of the default modehypothesis.
Proceedings of the National Academy of Sciences of theUnited States of America , 253–8 (2003). URL .[8] Galan, R. F. On How Network Architecture Determines the DominantPatterns of Spontaneous Neural Activity.
PLoS One , 1–10 (2008).[9] Honey, C. J., K¨otter, R., Breakspear, M. & Sporns, O. Network struc-ture of cerebral cortex shapes functional connectivity on multiple timescales. Proceedings of the National Academy of Sciences , 10240–10245 (2007). URL . .[10] Deco, G., Jirsa, V. K. & McIntosh, A. R. Emerging concepts for thedynamical organization of resting-state activity in the brain. NatureReviews Neuroscience , 43–56 (2011).1711] Deco, G., Senden, M. & Jirsa, V. How anatomy shapes dynamics :a semi-analytical study of the brain at rest by a simple spin model. Frontiers in Computation Neuroscience , 1–7 (2012).[12] Deco, G. et al. Resting-state functional connectivity emerges from struc-turally and dynamically shaped slow linear fluctuations.
The Journal ofNeuroscience , 11239–11252 (2013).[13] Goni, J. et al. Resting-brain functional connectivity predicted by ana-lytic measures of network communication.
Proceedings of the NationalAcademy of Sciences , 833–838 (2014). URL .[14] Sui, J. et al.
NeuroImage Three-way ( N-way ) fusion of brain imag-ing data based on mCCA + jICA and its application to discrimi-nating schizophrenia.
NeuroImage , 119–132 (2013). URL http://dx.doi.org/10.1016/j.neuroimage.2012.10.051 .[15] Sui, J. et al. Archival Report In Search of Multimodal NeuroimagingBiomarkers of Cognitive Deficits in Schizophrenia.
Biological Psychiatry , 794–804 (2015). URL http://dx.doi.org/10.1016/j.biopsych.2015.02.017 .[16] O’Muircheartaigh, J. & Jbabdi, S. Concurrent white matter bundlesand grey matter networks using independent component analysis. Neu-roImage , 296–306 (2018). URL https://doi.org/10.1016/j.neuroimage.2017.05.012 .[17] Calamante, F., Smith, R. E., Liang, X., Zalesky, A. & Connelly,A. Track-weighted dynamic functional connectivity (TW-dFC): a newmethod to study time-resolved functional connectivity.
Brain Structureand Function , 3761–3774 (2017).[18] Robinson, P. A. Interrelating anatomical, effective, and functional brainconnectivity using propagators and neural field theory.
Physical ReviewE , 1–18 (2012).[19] Robinson, P. A. et al.
Eigenmodes of brain activity: Neural field theorypredictions and comparison with experiment.
NeuroImage , 79–98(2016). 1820] Atasoy, S., Donnelly, I. & Pearson, J. Human brain networks functionin connectome-specific harmonic waves.
Nature communications , 1–10(2016).[21] Gu, S. et al. Controllability of structural brain networks.
Nature com-munications (2015).[22] Huang, W. et al. A Graph Signal Processing Perspective on FunctionalBrain Imaging.
Proceedings of the IEEE , 868–885 (2018).[23] Buades, A., Coll, B. & Morel, J.-M. A review of image denoising al-gorithms, with a new one.
SIAM Journal on Multiscale Modeling andSimulation: A SIAM Interdisciplinary Journal , 490–530 (2005).[24] Rudin, L. I., Osher, S. & Fatemi, E. Nonlinear total variation basednoise removal algorithms. Physica D: Nonlinear Phenomena , 259 –268 (1992).[25] Chambolle, A. An algorithm for total variation minimization and ap-plications. J. Math. Imaging Vis. , 89–97 (2004). URL http://dx.doi.org/10.1023/B:JMIV.0000011325.36760.1e .[26] Cands, E., Romberg, J. & Tao, T. Stable signal recovery from incompleteand inaccurate measurements. Communications on Pure and AppliedMathematics , 1207–1223 (2006). URL https://onlinelibrary.wiley.com/doi/abs/10.1002/cpa.20124 .[27] Shuman, D. I., Narang, S. K., Frossard, P., Ortega, A. & Vandergheynst,P. The emerging field of signal processing on graphs: Extending high-dimensional data analysis to networks and other irregular domains. IEEE Signal Processing Magazine , 83–98 (2013).[28] Chen, S., Sandryhaila, A., Moura, J. M. F. & Kovaˇcevi´c, J. SignalRecovery on Graphs. IEEE Transactions on Signal Processing , 1–20(2014). URL http://arxiv.org/abs/1411.7414 .[29] Koch, M. & Norris, D. An investigation of functional and anatomicalconnectivity using magnetic resonance imaging. Neuroimage , 241–250 (2002). 1930] Biswal, B., Yetkin, F., Haughton, V. & JS, H. Functional connectivity inthe motor cortex of resting human brain using echo-planar MRI. MagnReson Med , 537–541 (1995).[31] Liu, X. & Duyn, J. H. Time-varying functional network informationextracted from brief instances of spontaneous brain activity. Proceedingsof the National Academy of Sciences , 4392–4397 (2013). URL .[32] Andrews-Hannah, J. R., Reidler, J. S., Sepulcre, J., Poulin, R. & Buck-ner, R. L. Functional-anatomic fractionation of the brain’s default net-work.
Neuron , 550 – 562 (2010).[33] Zeithamova, D., Dominick, A. L. & Preston, A. R. Article Hippocam-pal and Ventral Medial Prefrontal Activation during Retrieval-MediatedLearning Supports Novel Inference. Neuron
Ann NY Acad Sci. , 1–38 (2008).[35] Fox, M. D. et al. The human brain is intrinsically organized into dy-namic, anticorrelated functional networks.
Proceedings of the NationalAcademy of Sciences , 9673 – 9678 (2005).[36] Greicius, M. D., Srivastava, G., Reiss, A. L. & Menon, V. Default-modenetwork activity distinguishes Alzheimer’ s disease from healthy aging :Evidence from functional MRI.
Proceedings of the National Academy ofSciences , 4637–4642 (2004).[37] Piccoli, T., Valente, G., Linden, D. E. J. & Re, M. The Default ModeNetwork and the Working Memory Network Are Not Anti-Correlatedduring All Phases of a Working Memory Task.
PLOS One
Nature Communications (2015).[39] Leech, R., Kamourieh, S., Beckmann, C. F. & Sharp, D. J. Fractionatingthe Default Mode Network : Distinct Contributions of the Ventral and20orsal Posterior Cingulate Cortex to Cognitive Control. The Journal ofNeuroscience , 3217–3224 (2011).[40] Squire, L., Genzel, L., Wixted, J. & Morris, R. Memory Consolidation. Cold Spring Harbor Perspectives in Biology , 1–21 (2015).[41] Hanlon, F. M. et al. Frontotemporal anatomical connectivity andworking-relational memory performance predict everyday functioning inschizophrenia.
Psychophysiology , 1340–1352 (2012).[42] Wheeler-kingshott, C. A. M. & Cercignani, M. About Axial and RadialDiffusivities. Magn Reson Med , 1255–1260 (2009).[43] Maier-Hein, K. H. et al.
The challenge of mapping the human connec-tome based on diffusion tractography.
Nature communications , 1349(2017).[44] Behjat, H., Leonardi, N., S¨ornmo, L. & Van De Ville, D. Anatomically-adapted graph wavelets for improved group-level fMRI activation map-ping. Neuroimage , 185–199 (2015).[45] Logothetis, N. & Wandell, B. Interpreting the BOLD signal.
Annu. Rev.Physiol. , 735–769 (2004).[46] Gawryluk, J. R., Mazerolle, E. L., Arcy, R. C. N. D. & Chan, K. C. Doesfunctional MRI detect activation in white matter ? A review of emergingevidence , issues , and future directions. Frontiers in Neuroscience ,1–12 (2014).[47] Peer, X. M., Nitzan, M., Bick, A. S., Levin, N. & Arzy, S. Evidencefor Functional Networks within the Human Brain ’ s. Journal of Neu-roscience , 6394–6407 (2017).[48] Ding, Z. et al. Detection of synchronous brain activity in white mat-ter tracts at rest and under functional loading.
Proceedings of the Na-tional Academy of Sciences of the United States of America , 595–600 (2018).[49] Fraser, L. M., Stevens, M. T., Beyea, S. D. & Arcy, R. C. N. D. Whiteversus gray matter : fMRI hemodynamic responses show similar char-acteristics , but differ in peak amplitude.
BMC Neuroscience (2012).2150] Newton, A. T., Anderson, A. W., Ding, Z. & Gore, J. C. Character-ization of the hemodynamic response function in white matter tractsfor event-related fMRI.
Nature Communications http://dx.doi.org/10.1038/s41467-019-09076-2 .[51] Beckmann, C., DeLuca, M., Devlin, J. & Smith, S. Investigations intoresting-state connectivity using independent component analysis.
PhilosTrans R Soc Lond B Biol Sci , 1001–1013 (2005).[52] Leonardi, N. et al.
NeuroImage Principal components of functional con-nectivity : A new approach to study dynamic brain connectivity duringrest.
NeuroImage , 937–950 (2013). URL http://dx.doi.org/10.1016/j.neuroimage.2013.07.019 .[53] Markov, N. T., Lamy, C., Essen, D. C. V., Knoblauch, K. & Lyon, D. APredictive Network Model of Cerebral Cortical Connectivity Based ona Distance Rule. Neuron , 184–197 (2012).[54] Alexander-bloch, A. F. et al. The Anatomical Distance of Func-tional Connections Predicts Brain Network Topology in Health andSchizophrenia.
Cerebral Cortex et al.
Functional alignment with anatomical networks isassociated with cognitive flexibility.
Nature Human Behaviour , 156–164 (2018). URL https://doi.org/10.1038/s41562-017-0260-9 .[56] Park, H.-j. et al. Corpus Callosal Connection Mapping Using CorticalGray Matter Parcellation and DT.
Human Brain Mapping , 503–516(2008).[57] Glasser, M. F. et al.
NeuroImage The minimal preprocessing pipelinesfor the Human Connectome Project.
NeuroImage , 105–124 (2013).[58] Van Dijk, K. R. et al. Intrinsic functional connectivity as a tool forhuman connectomics: theory, properties, and optimization.
Journal ofneurophysiology , 297–321 (2009).22 cknowledgements
We thank T. Bolton for his assistance and valuable comments. This workwas supported by the Swiss National Science Foundation under the ProjectGrant 205321-163376.
Author Contributions
A.T., D.A., and H.B. designed the construction of the voxel-wise braingrid, A.T. and D.V.D.V. designed and performed research. A.T and D.V.D.V.wrote the manuscript.
Competing Interests
The authors declare that they have no competing financial interests.