An algebraic topological method for multimodal brain networks comparisons

Understanding brain connectivity is one of the most important issues in neuroscience. Nonetheless, connectivity data can reflect either functional relationships of brain activities or anatomical connections between brain areas. Although both representations should be related, this relationship is not straightforward. We have devised a powerful method that allows different operations between networks that share the same set of nodes, by embedding them in a common metric space, enforcing transitivity to the graph topology. Here, we apply this method to construct an aggregated network from a set of functional graphs, each one from a different subject. Once this aggregated functional network is constructed, we use again our method to compare it with the structural connectivity to identify particular brain regions that differ in both modalities (anatomical and functional). Remarkably, these brain regions include functional areas that form part of the classical resting state networks. We conclude that our method -based on the comparison of the aggregated functional network- reveals some emerging features that could not be observed when the comparison is performed with the classical averaged functional network.

1. Introduction

In the last decade, the use of advanced tools derived from neuroimaging and complex networks theory have significantly improved our understanding of brain functioning (Sporns, 2011). Notably, connectivity-based methods have had a prominent role in characterizing normal brain organization as well as alterations due to various brain disorders (Varela et al., 2001; Stam and van Straaten, 2012; Stam, 2014). Most of the recent works aim to quantify the role of connectivity in the communication abilities of neural systems. However, the very same notion of connectivity is controversial since data used in brain connectivity studies can reflect functional neural activities (electrical, magnetic or hemodynamic/metabolic) or anatomical properties (Varela et al., 2001; Bullmore and Sporns, 2009). Neuroanatomical connectivity is meant as the description of the physical connections (axonal projections) between two brain sites (Bullmore and Sporns, 2009), whereas functional connectivity is defined as the estimated temporal correlation between spatially distant neurophysiological activities such as electroencephalographic (EEG), magnetoencephalographic (MEG), functional magnetic resonance imaging (fMRI) or positron emission tomography (PET) recordings (Varela et al., 2001).

In recent years, the concept of “brain networks” is becoming fundamental in neuroscience (Stam and Reijneveld, 2007; Bullmore and Sporns, 2009; Stam and van Straaten, 2012; Stam, 2014). Within this framework, nodes stand for different brain regions (e.g., parcelated areas or recording sites) and links indicate either the presence of an anatomical path between those regions or a functional dependence between their activities. In the last years, this representation of the brain has allowed to visualize and describe its non-trivial topological properties in a compact and objective way. Nowadays, the use of network-based analysis in neuroscience has become essential to quantify brain disfunctions in terms of aberrant reconfiguration of functional brain networks (Stam and Reijneveld, 2007; Stam and van Straaten, 2012; Stam, 2014).

Experimental evidence has revealed, for instance, alterations in functional and anatomical brain networks in normal cognitive processes, across development, and in a wide range of neurological diseases (see Bullmore and Sporns, 2009; Stam, 2014; and references therein). Despite its evident interplay, comparison between anatomical and functional brain networks is not straightforward (Deco et al., 2011; Nicosia et al., 2014). Theoretical studies provide support for the idea that structural networks determine some aspects of functional networks (Deco et al., 2011), but it is less clear how the anatomical connectivity supports or facilitates the emergence of functional networks. Although nodes with similar connection patterns tend to exhibit similar functionality, the functionality of an individual neural node is strongly determined by the pattern of its interconnections with the rest of the network (Nicosia et al., 2014).

Correspondences between functional and structural networks remains thus an active research area (Honey et al., 2007, 2009, 2010). A better understanding of how anatomical scaffolds support functional communication of brain activities is necessary to better understand normal neural processes, as well as to improve identification and prediction of alterations in brain diseases.

In this paper we address this relationship between anatomical and functional connectivity. In previous studies, the correspondence of these networks has been often assessed by the difference in an Euclidean space of vectors containing connectivity measures such as the clustering coefficient, shortest path length, degree distribution, etc. Here, we propose a radically different framework for studying brain connectivity differences. Instead of extracting a vector of features for each network (anatomical or functional), we jointly embed all of them in a common metric space that allow straightforward comparisons. Before embedding the functional and the anatomical networks into the common metric space, we aggregate a group of subjects (e.g., functional networks) according to Simas et al., (submitted) to obtain a group representation network. The method employed in this work allows to preserve connected components and to identify, among different subjects, a common underlying network structure. Our approach may provide a useful insight for the analysis of multiple networks obtained from multiple brain modalities or groups (healthy volunteers vs. patients, for instance).

2. Methods and Materials

2.1. fMRI and DTI Data

In this study we consider anatomical and functional brain connectivities-extracted from diffusion-weighted DW-MRI and fMRI data, respectively- defined on the same brain regions. Brain images were partitioned into the 90 anatomical regions (N = 90 nodes of the networks) of the Tzourio-Mazoyer brain atlas (Tzourio-Mazoyer et al., 2002) using the automated anatomical labeling method.

The anatomical connectivity network is based on the connectivity matrix obtained by Diffusion Magnetic Resonance Imaging (DW-MRI) data from 20 healthy participants, as described in Iturria-Medina et al. (2008). The elements of this matrix represent the probabilities of connection between the 90 brain regions of interest. These probabilities are proportional to the density of axonal fibers between different areas, so each element of the matrix represents an approximation of the connection strength between the corresponding pair of brain regions.

The functional brain connectivity was extracted from BOLD fMRI resting state recordings obtained as described in Valencia et al. (2009). All acquired brain volumes were corrected for motion and differences in slice acquisition times using the SPM5¹ software package. All fMRI data sets (segments of 5 min recorded from healthy subjects) were co-registered to the anatomical data set and normalized to the standard MNI (Montreal Neurological Institute) template image, to allow comparisons between subjects. As for DW-MRI data, normalized and corrected functional scans were sub-sampled to the anatomical labeled template of the human brain (Tzourio-Mazoyer et al., 2002). Regional time series were estimated for each individual by averaging the fMRI time series over all voxels in each of the 90 regions. To eliminate low frequency noise (e.g., slow scanner drifts) and higher frequency artifacts from cardiac and respiratory oscillations, time-series were digitally filtered with a finite impulse response (FIR) filter with zero-phase distortion (bandwidth 0.01—0.1 Hz) as in Valencia et al. (2009).

A functional link between two time series x_i(t) and x_j(t) (normalized to zero mean and unit variance) was defined by means of the linear cross-correlation coefficient computed as r_ij = 〈x_i(t)x_j(t)〉, where 〈·〉 denotes the temporal average. For the sake of simplicity, we only considered here correlations at zero lag. To determine the probability that correlation values are significantly higher than what is expected from independent time series, r_ij(0) values (denoted r_ij) were firstly variance-stabilized by applying the Fisher's Z transform.

$\begin{array}{ll}{Z}_{ij}=0.5\ln \left(\frac{1+r_{ij}}{1-r_{ij}}\right)\hfill & \left(1\right)\hfill \end{array}$

Under the hypothesis of independence, Z_ij has a normal distribution with expected value 0 and variance 1/(df_ij−3), where df is the effective number of degrees of freedom (Bartlett, 1946; Bayley and Hammersley, 1946; Jenkins and Watts, 1968). If the time series consist of independent measurements, df_ij simply equals the sample size, N. Nevertheless, autocorrelated time series do not meet the assumption of independence required by the standard significance test, yielding a greater Type I error (Bartlett, 1946; Bayley and Hammersley, 1946; Jenkins and Watts, 1968). In presence of auto-correlated time series df must be corrected by the following approximation $\frac{1}{d{f}_{ij}}\approx \frac{1}{N}+\frac{2}{N}{\displaystyle {\sum }_{\tau }\frac{N-\tau }{N}{r}_{ii}(\tau )r_{jj}(\tau )}$, where r_xx(τ) is the autocorrelation of signal x at lag τ.

2.2. Networks Normalization

From Equation (1) our networks weights are in a non-normalized interval Z_ij ∈ [a,b] ⊂ ℝ. In order to apply the framework described in Simas and Rocha (2015), we normalize our networks weights into the unit interval I = [0, 1] by means of a unique linear function:

$\begin{array}{ll}{w}_{ij}=\frac{(1-2ϵ)Z_{ij}+(2ϵ-1)·MIN(Z_{ij})}{MAX(Z_{ij})-MIN(Z_{ij})}+ϵ\hfill & \left(2\right)\hfill \end{array}$

where ϵ in general is set to 0.01 in order to avoid merging and isolate vertices with weights at the boundaries of Z_ij ∈ [a,b]. As proved in Simas and Rocha (2015), since the normalization is done by a unique linear function this does not affect networks properties.

2.3. fMRI Networks Aggregation and Embedding

Among many ways to aggregate a group of networks here we employ a topological algebraic way to aggregate a group of networks. The networks group possess the same nodes but different edges values and can mathematical be represented by a weighted graph G = (N,E). N is the set of nodes representing the brain ROI's (N = 90 in this study) and E is the set of edges values (connections) between ROI's, e.g., ∀e_i,j ∈ E: e_i,j ∈ [0, 1] in the proximity space or ∀d_i,j ∈ E: d_i,j ∈ ℝ⁺₀ ∪ {0, + ∞} in the distance space.

For the sake of simplicity, we denote a network with the same notation we use for the set of nodes N, i.e., a set of n networks (e.g., group of subjects) is represented by {N_k} with k∈ {1, 2, 3.…, n}.

One possible way to aggregate a group of n networks is simply by averaging the homologous edges values. Obtaining in this way a group representative network, N*.

$\begin{array}{ll}{N}_{i,j}^{\ast }=e_{i,j}^{\ast }=\frac{{\displaystyle {\sum }_{k=1}^n{e}_{i,j}^{[k]}}}{n}\hfill & \left(3\right)\hfill \end{array}$

where e^[k]_i,j is the edge e_i,j from network N_k.

Another way to aggregate networks, as explained in Simas et al., (submitted), is by considering all networks as a multilayer network (often called multiplex), which can be represented as a fourth-order tensor Simas et al., (submitted). This tensor can be represented as a extended matrix (Sole-Ribalta et al., 2013). The work of Simas and Rocha (2015), introduces a framework to aggregate networks in an algebraic way, relating it with fuzzy logic reasoning, and in Simas et al., (submitted) this work was extended to multilayer networks. In order to work algebraically with networks we have to set an algebra (defined as a vector space equipped with a bilinear product). This algebra allows us to perform algebraic operations with networks in the same way we perform algebraic operations in other contexts with other algebras (such as adding and multiplying real numbers). In short, a network can be represented by an adjacency matrix and a multilayer network by a tensor. Considering a set of tensors working under the algebra L = (I, ⊕, ⊗), where the weights (tensor entries) of the tensors in I ⊆ ℝ (subset of extended real line) and ⊕ and ⊗ two binary operators, we can represent a multilayer network with tensor T in this algebra. In Simas et al., (submitted) we have shown the particular case of multiplex networks, where layers are connected with weights w_i,i,L_k,L_j = 1 (in the proximity space), that the representative group network (e.g., functional) can be represented by N* in the distance space (see below and Equation 6), as:

$\begin{array}{ll}{N}^{\ast }\equiv N_1\oplus N_2\oplus ···\oplus N_k\hfill & \left(4\right)\hfill \end{array}$

and the respective embedding by the following equation:

$\begin{array}{ll}{N}_{embedded}=N^{\ast }\oplus N^{\ast 2}\oplus ···\oplus N^{\ast r}\hfill & \left(5\right)\hfill \end{array}$

where N* is defined in Equation (4) and r, is the convergence parameter (Simas and Rocha, 2015; Simas et al., submitted). Figure 1 summarizes the metric embedding of a multiplex network described above.

FIGURE 1

Figure 1. Schematic representation of the main steps for the described networks aggregation and metric embedding (defined here for the algebra L).

Embedding a network of networks or, in our specific case, a multiplex fMRI network, allows us to determine which edges in the several layers contribute to the aggregation. We can therefore determine the subjects that contribute more/less or none to the aggregated network, and identify in each subject the sub-graphs for which they may have the highest contribution.

For our particular case, we embed our networks using the Metric Closure (Simas and Rocha, 2015) defined by the algebra L = (ℝ⁺₀ ∪ {+∞}, min, +), where ⊕ = min and ⊗ = +. The metric closure or metric embedding of a given network into a metric space, is a generalization of All Pairs Shortest Paths Problem (APSP) as shown in Simas and Rocha (2015). In this case the Johnson algorithm can be used to calculate the metric closure (Johnson, 1977).

Note that to calculate the metric closure (based on the Johnson algorithm) of a network we have to translate our networks from a proximity space into a distance space. There are many possible mappings to map a similarity space into a distance space, see (Simas and Rocha, 2015). Applying Equation (6) to all network weights, w_i,j ∈ [0, 1] (for more details see Simas and Rocha, 2015), we obtain the isomorphic distance network with weights d_i,j ∈ ℝ⁺₀ ∪ {+∞}.

$\begin{array}{ll}{d}_{i,j}=\frac{1}{w_{i,j}}-1\hfill & \left(6\right)\hfill \end{array}$

The formalism behind the metric closure should not be confused with the formalism in Tropical Algebra geometry (Pachter and Sturmfels, 2004; Theobald, 2006). Both formalisms employ the same algebra for the isomorphism d_i,j = φ(x) = −log(x), which corresponds to a Schweizer-Sklar or Frank t-norm generator with λ = 0 or 1 (see Equation 6 and Simas and Rocha, 2015) under the formalism in Simas and Rocha (2015). The formalism in Simas and Rocha (2015) uses any isomorphism φ to set a specific metric into a weighted graph when translated to the isomorphic distance space. A more detailed discussion on this relation between the work (Simas and Rocha, 2015) and Tropical geometry can be found in Simas et al., (submitted).

Embedding networks or multilayer networks allows us: (a) to detected clusters of nodes in a high-dimensional topological spaces, and by projecting the algebraic high-dimensional embedding into 3D, (b) it allows to perform exploratory networks analysis (c) to preserve the multilayer sub-structures across layers/subjects, better than other aggregations methods, as compared with the specific case of “simple” averaging (Equation 3).

Next we compare both methods of aggregation, “simple” averaging (Equation 5) and algebraic aggregation (Equation 5 according to Simas et al., submitted) of our fMRI networks, respectively, using our proposed method of embedding and comparing networks.

2.4. Multimodal Networks Comparison

In general, networks have been compared using statistical measures of local and global properties, such as: clustering coefficient, small-worldness, degree distributions, etc. We can find in the literature some examples of such techniques to compare multiple networks (Bullmore and Sporns, 2009; Stam, 2014). Our approach in this work is different. After embedding networks into the same metric space defined by the applied algebra, in our case L = (ℝ⁺₀ ∪ {+ ∞}, min, +), we are able to compared them topologically. However, since networks generally come from different modalities (e.g., fMRI and DTI) it requires a previous step. We need to normalize the embedded edge weights distributions from the different modalities to the same average and variance to remove scale factors. One possible way to normalize both distributions, if we assume normality, is by calculating the z-score of the edge weights distributions (zero average and standard deviation set to the unit).

The embedded networks represent a hyper-grid in a multi-dimensional space with dimension equal or below to the number of nodes. In order to simplify and have some visual insight we can downgrade linearly this multidimensional grid into a 3D grid. This can be achieved applying to the embedded networks any technique for dimensionality reduction such as linear/non-linear Multi-Dimensional Scaling (MDS). MDS procedures refer to a set of related ordination techniques used in information visualization, in particular to display the information contained in a distance matrix (Borg and Groenen, 2005). These techniques guarantee, with a given distortion, that the relative distance between nodes is preserved in both multi-dimensional and low-dimensional reduction space. Plotting this low-dimensional grid (e.g., in 3D) we can use any statistical technique to fit a continuous surface into the data (see below Figures 3, 4). Its is natural to think that the difference between two surfaces obtained from different networks will emphasize topologically differences between the two connectivities. In this work we performed this operation in the multi-dimensional space by subtracting homologous embedded edges weights and take the absolute value of both embedded hyper-grids. This is, we subtract homologous embedded edges pairwise according to the formula:

$\begin{array}{ll}\begin{array}{l}M=|M_{fMR{I}^{\ast }}-M_{DTI}|\equiv \\ \{e_{i,j}^{diff}:\forall e_{i,j}^{\ast }\in M_{fMR{I}^{\ast }}\wedge \forall e_{i,j}^{\ast \ast }\in M_{DTI}:e_{i,j}^{diff}=|e_{i,j}^{\ast }-e_{i,j}^{\ast \ast }\left|\right\}\end{array}\hfill & \left(7\right)\hfill \end{array}$

M is the difference grid in the multi-dimensional space. Because the M-grid represents the difference between the two grids from different modalities (see above), the relative distance between nodes in M (given by Equation 7) should be concentrated at the origin if they are topological similar, otherwise widely distributed in the multi-dimensional space. Nodes at a distance from the origin of s-standard deviations are statistical different. Moreover, since we z-scored both embedded edge distributions this give us some degree of statistical significance when we compare both networks. All nodes that lay outside of a hyper-sphere with center at the origin with radius R = s, are statistical different. Here we had set σ = 1 for both distributions (z-score variables are estimated from the distributions of the embedded weights).

Figure 2 illustrates this process. After applying to both fMRI* and DTI networks the same algebra and the metric embedding described above, both networks rely on the same metric space, therefore comparable. Topological differences can be visually seen in a linearly downgraded to 3D dimensions using a MDS technique, which preserves the relative distance between points in the grid (nodes or brain areas).

FIGURE 2

Figure 2. Topological algebraic networks comparison. Connectivity from different modalities (here fMRI and DTI) are firstly embedded (black dot points on the manifolds indicate the brain nodes) and then compared in a low-dimensional space. Black points outside the sphere correspond to nodes with a topological difference (at a given threshold) in the two modalities.

3. Results

In Figure 3 we illustrate the results of different aggregation procedures on the ensemble of fMRI networks. Compared with a fMRI connectivity matrix from a single subject (Figure 3A), one can notice the difference of a single averaging across subjects (Figure 3B) and our proposed algebraic topologically aggregated connectivity network (Figure 3C). It is clear that the averaging procedure tends to blur connectivity values between nodes. In contrast, the topologically algebraic aggregation can preserve components that are common across subjects. As other multilinear algebra or tensor-based analysis, our approach provide a natural mathematical framework for studying connectivity data with multidimensional structure. For illustrative purposes, we also show the DTI connectivity matrix in Figure 3D). It worths noticing the similarity of the anatomical connectivity structure with the aggregated (multiplex) connectivity obtained in Figure 3C. Moreover, since each layer encodes the functional network for a given subject, each subject contributes to the tensor aggregation/embedding with some or none connections (edges), as depicted in metric closure, Equation (5). If a layer do not contribute for the aggregation/embedding, we may consider this layer (subject network) as an outlier. Moreover, we are also able to identify the specific sub-network contribution (edges) of a given layer to the aggregation/embedding.

FIGURE 3

Figure 3. (A) fMRI single subject network (B) Average aggregated fMRI network (C) fMRI Algebraic Topologically aggregated (multiplex) network (D) DTI network.

Low-dimensional embeddings of different aggregated networks are illustrated in Figure 4. High-dimensional data, such as the information contained in the distance matrix obtained for the different networks, can be difficult to interpret. Here, multidimensional scaling (MDS) was used for visualizing the level of similarity of individual nodes of each -aggregated- network. The MDS algorithm aims to place each node in a low dimensional space such that the between-nodes distances are preserved as much as possible. This representation into a low-dimensional space enables an exploratory analysis and makes data analysis algorithms more efficient. Indeed, from the different plots of Figure 4 one can identify brain areas that are topologically close in the aggregated network as those points that are close on the 3D grid. This is clearly illustrated by the MDS representation of the multiplex functional network (Figure 4C). Nodes from the occipital regions form a compact group of nodes topologically close (with similar connectivity structure), as revealed by the blue points depicted on Figure 4C. We also notice that a compact group of nodes is formed by regions of the temporal lobe, putamen and insula, which are indicated by the red circle. Similarly, the anatomical network in Figure 4D clearly displays a natural organization, i.e., nodes of the two hemispheres lie on both sides of the dotted black line. Further, nodes from occipital regions in the anatomical network, indicated by the blue circles (including calcarina, cuneus, precuneus, …), are distantly located from the group of frontal brain areas indicated by the red marks.

FIGURE 4

Figure 4. Multi-Dimentional Scaling (MDS) of the embedded networks (A) fMRI single subject (B) fMRI average embedded network (C) fMRI Algebraic Topological aggregation (multiplex) embedded network (D) DTI embedded network. Black dots indicate the embedded nodes. In plots (C,D), blue and red points indicate the groups of brain areas discussed in the text.

Finally, Figure 5 displays the difference grid M in a low-dimensional space. As defined in Equation (7), M corresponds to the relative distance between nodes in networks from different modalities. Differences between brain areas are represented as points widely distributed in the low-dimensional space. Those nodes from different modalities (fMRI and DTI) that share an identical topological structure are located at the origin. The larger the difference in the connectivity structure, the larger the distance from the origin. By setting a threshold s, one can identify brain areas with similar connectivity as those points that lie inside of the hyper-sphere of radius s with center at the origin.

FIGURE 5

Figure 5. Comparisons between DTI and all other embedded fMRI networks. (A) 3D projections from Equation (7). Only points outside the sphere are plotted. (B) Number of ROI's inside the sphere of radius of s. Results from a single subject, average connectivity and multiplex networks are represented by the red, blue, and green points and curves, respectively. We consider the regions statistically different for s>1 and statistically equal for s<1. This shows that the multiplex algebraic aggregation (green) is more similar algebraically to DTI then average aggregation (blue) and single subject fMRI network (red).

The number of brain regions (ROIs) with similar anatomical and functional connectivity are given in Figure 5B as a function of the threshold s. Curves correspond to the number of regions inside a hyper-sphere of various radius. We notice that the number of regions differ as a function of the aggregated network's type. It is worthy to mention that the differences above s-standard deviations are the important ones, since is above this threshold that the ROI's or nodes become statistical different when compare networks. In our example, the fluctuations below one standard deviations may give us some trend but all nodes in the networks are statistical equal for all types of aggregation. For our specific case, as an example, the brain areas located outside the hyper-sphere of radius s = 1.2 for the two types of aggregation, are listed in the Table 1.

TABLE 1

Table 1. ROIs with connectivity differences from DTI at 1.2 standard deviation.

4. Discussion

The recent prevalence of applications involving multidimensional and multimodal brain data has increased the demand for technical developments in the analysis of such complex data. Indeed, the discrepancy between structural and functional brain connectivity is a current challenge for understanding general brain functioning. In this paper, we presented a method for characterizing the correspondences between functional and anatomical connectivity. To summarize, the main steps of our method are:

1. Metric network embedding: This procedure embed a group of connectivity graphs in a common space allowing straightforward comparisons. In contrast with simple averaging of connectivity matrices, the topologically algebraic aggregation can preserve components that are common across different subjects or different neuroimaging modalities. This tensor-based aggregation allows enhancing the common underlying structures providing a natural mathematical framework for studying connectivity data with multidimensional structure.

2. Multimodal Networks comparison: the differences between the embedded networks are calculated and represented in a low-dimensional space. Multi-Dimensional Scaling simply enables to display the information contained in the resulting distance matrix allowing thus an exploratory analysis of the data.

3. Detection of nodes (ROI's) with different connectivities: from points widely distributed in the low-dimensional space one can detect brain nodes that share a similar topological structure as those points are located close to the the origin. One can identify brain areas with the largest difference between anatomical and functional connectivity as those points located outside an imaginary hyper-sphere of a radius given by a threshold (Table 1)

Our findings suggest that embedding a brain network on a metric space may reveal regions that are members of large areas or subsystems rather than regions with a specific role in information processing. This is clearly illustrated for the anatomical network in Figure 5D, where frontal and occipital brain areas of both hemispheres are situated at distantly and located points of the space. Contrary to a classical averaging of connectivity matrices, the embedding of the multiplex functional network reveals brain areas that play a role in large brain system such as the occipital regions, known to be active when the subject is at wakeful rest.

Although experimental evidence suggests that functionally linked brain regions have an underlying structural core, this relationship does not exhibit a simple one-to-one mapping (Wang et al., 2014). These correspondences have also been investigated in specific subsystems, must of them focused on the default mode network (DMN), which is a group of brain regions that preferentially activate when individuals engage in internal tasks, i.e., when the subject is not focused on the outside world but the brain is at wakeful rest. Several studies report that the DMN exhibits a high overlap in its structural and functional connectivity (Honey et al., 2009; Wang et al., 2014). Nevertheless, strong discrepancies have been reported and strong functional links can be found between regions without direct structural linkages (Honey et al., 2009).

At a group level, one of the reasons for this discrepancy between structural and functional connectivity has been suggested to be the functional variability across subjects (Skudlarski et al., 2008; Honey et al., 2009; Wang et al., 2014). Indeed, clinical studies have provided evidence for a large heterogeneity of the functional connectivity, particularly in groups of patients with brain disorders such as neuropsychiatric disorders, which strongly alters the structural-functional relationships (Wang et al., 2014). Analytical tools are therefore required to account for this variability in order to enhance the common underlying network structure.

Results suggest that averaged aggregation captures the general differences in regions that play a role in visual, auditory and body self-awareness processes, but fails to identify in detail other specific areas across the subjects/groups. In Table 1 we observed that the average aggregation essentially captures part of visual (calcarine, cuneus, lingual, occipital), auditory (superior temporal gyrus), and insula regions that are associated to visual process and body self-awareness. Detection of visual and auditory regions suggest that the averaged aggregation mainly capture regions activated by the resting state condition of the recordings.

From the multiplex aggregation (or algebraic aggregation) shown in Table 1, we observed that besides capturing the well-known visual (occipital areas), primary sensory cortex (postcentral), and auditory regions (Heschl gyrus, superior temporal, thalamus), this approach also captures some other network sub-structures involved in touch activation (postcentral gyrus, thalamus) and emotional state activations (amygdala, thalamus, posterior cingulate). This alines with our claim that algebraic aggregation preserves better the multilayer sub-structures across a group of subjects (multilayers) accounting for as much of the variability in the data as possible.

Although we cannot definitively provide a one-to-one mapping of the structural and functional connectivity, we think that our method could provide new insights on the organization of brain networks during diverse cognitive or pathological states. We therefore hope that our approach will foster more principled and successful analysis of multimodal brain connectivity datasets.

For all the methods described in this article we provide the corresponding MATLAB software code. Data and code are freely available at the website https://sites.google.com/site/fr2eborn/download.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors thank Y. Iturria-Medina and L. Melie-García for sharing the DTI connectivity data used in the study. This work was supported by the LASAGNE (Contract No.318132) and MULTIPLEX (Contract No.317532) EU projects. AD acknowledges support from Generalitat de Catalunya (2014SGR608) and Spanish MINECO (FIS2012-38266).

Footnotes

1. ^http://www.fil.ion.ucl.ac.uk

References

Bartlett, M. S. (1946). On the theoretical specification and sampling properties of autocorrelated time-series. Supp. J. R. Stat. Soc. 8, 27–41. doi: 10.2307/2983611

CrossRef Full Text | Google Scholar

Bayley, G. V., and Hammersley, J. M. (1946). The effective number of independent observations in an autocorrelated time series. Supp. J. R. Stat. Soc. 8, 184–197. doi: 10.2307/2983560

CrossRef Full Text | Google Scholar

Borg, I., and Groenen, P. (2005). Modern Multidimensional Scaling: Theory and Applications. New York, NY: Springer-Verlag.

Google Scholar

Bullmore, E., and Sporns, O. (2009). Complex brain networks: graph theoretical analysis of structural and functional systems. Nat. Rev. Neurosci. 10, 186–198. doi: 10.1038/nrn2575

PubMed Abstract | CrossRef Full Text | Google Scholar

Deco, G., Jirsa, V. K., and McIntosh, A. R. (2011). Emerging concepts for the dynamical organization of resting-state activity in the brain. Nat. Rev. Neurosci. 12, 43–56. doi: 10.1038/nrn2961

PubMed Abstract | CrossRef Full Text | Google Scholar

Honey, C. J., Kötter, R., Breakspear, M., and Sporns, O. (2007). Network structure of cerebral cortex shapes functional connectivity on multiple time scales. Proc. Natl. Acad. Sci. U.S.A. 104, 10240–10245. doi: 10.1073/pnas.0701519104

PubMed Abstract | CrossRef Full Text | Google Scholar

Honey, C. J., Sporns, O., Cammoun, L., Gigandet, X., Thiran, J. P., Meuli, R., et al. (2009). Predicting human resting-state functional connectivity from structural connectivity. Proc. Natl. Acad. Sci. U.S.A. 106, 2035–2040. doi: 10.1073/pnas.0811168106

PubMed Abstract | CrossRef Full Text | Google Scholar

Honey, C. J., Thivierge, J. P., and Sporns, O. (2010). Can structure predict function in the human brain? Computational models of the brain. Euroimage 52, 766–776. doi: 10.1016/j.neuroimage.2010.01.071

PubMed Abstract | CrossRef Full Text | Google Scholar

Iturria-Medina, Y., Sotero, R. C., Canales-Rodríguez, E. J., Alemán-Gómez, Y., and Melie-García, L. (2008). Studying the human brain anatomical network via diffusion-weighted MRI and graph theory. Neuroimage 40, 1064. doi: 10.1016/j.neuroimage.2007.10.060

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenkins, G. M., and Watts, D. G. (1968). Spectral Analysis and Its Applications. San Francisco, CA: Holden-Day.

Johnson, D. B. (1977). Efficient algorithms for shortest paths in sparse networks. J. ACM 24, 1–13.

Google Scholar

Nicosia, V., Valencia, M., Chavez, M., Díaz-Guilera, A., and Latora, V. (2014). Remote synchronization reveals network symmetries and functional modules. Phys. Rev. Lett. 110:174102. doi: 10.1103/PhysRevLett.110.174102

PubMed Abstract | CrossRef Full Text | Google Scholar

Pachter, L., and Sturmfels, B. (2004). Tropical geometry of statistical models. Proc. Natl. Acad. Sci. U.S.A. 101, 16132–16137. doi: 10.1073/pnas.0406010101

PubMed Abstract | CrossRef Full Text | Google Scholar

Simas, T., and Rocha, L. M. (2015). Distance closures on complex networks. Netw. Sci. 3, 227–268. doi: 10.1017/nws.2015.11

CrossRef Full Text | Google Scholar

Skudlarski, P., Jagannathan, K., Calhoun, V. D., Hampson, M., Skudlarska, B. A., and Pearlson, G. (2008). Measuring brain connectivity: diffusion tensor imaging validates resting state temporal correlations. Neuroimage 43, 554–561. doi: 10.1016/j.neuroimage.2008.07.063

PubMed Abstract | CrossRef Full Text | Google Scholar

Sole-Ribalta, A., De Domenico, M., Kouvaris, N. E., Diaz-Guilera, A., Gomez, S., and Arenas, A. (2013). Spectral properties of the Laplacian of multiplex networks Phys. Rev. E 88:032807. doi: 10.1103/PhysRevE.88.032807

PubMed Abstract | CrossRef Full Text | Google Scholar

Sporns, O. (2011). Networks of the Brain. Boston, MA: Massachusetts Institute of Technology.

Google Scholar

Stam, C. J. (2014). Modern network science of neurological disorders. Nat. Rev. Neurosci. 15, 683–695. doi: 10.1038/nrn3801

PubMed Abstract | CrossRef Full Text | Google Scholar

Stam, C. J., and Reijneveld, J. C. (2007). Graph theoretical analysis of complex networks in the brain. Nonlin. Biomed. Phys. 1, 1753–4631. doi: 10.1186/1753-4631-1-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Stam, C. J., and van Straaten, E. C. W. (2012). The organization of physiological brain networks. Clin. Neurophysiol. 123, 1067–1087. doi: 10.1016/j.clinph.2012.01.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Theobald, T. (2006). On the frontiers of polynomial computations in tropical geometry. J. Symb. Comput. 41, 1360–1375. doi: 10.1016/j.jsc.2005.11.006

CrossRef Full Text | Google Scholar

Tzourio-Mazoyer, N., Landeau, B., Papathanassiou, D., Crivello, F., Etard, O., Delcroix, N., et al. (2002). Automated anatomical labeling of activations in SPM using a macroscopic anatomical parcellation of the MNI MRI single-subject brain. Neuroimage 15, 273–289. doi: 10.1006/nimg.2001.0978

PubMed Abstract | CrossRef Full Text | Google Scholar

Valencia, M., Pastor, M. A., Fernández-Seara, M. A., Artieda, J., Martinerie, J., and Chavez., M. (2009). Complex modular structure of large-scale brain networks. Chaos 19, 023119. doi: 10.1063/1.3129783

PubMed Abstract | CrossRef Full Text | Google Scholar

Varela, F., Lachaux, J.-P., Rodriguez, E., and Martinerie, J. (2001). The brainweb: phase synchronization and large-scale integration. Nat. Rev. Neurosci. 2, 229–239. doi: 10.1038/35067550

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Z., Dai, Z., Gong, G., Zhou, C., and He, Y. (2014). Understanding structural-functional relationships in the human brain: a large-scale network perspective. Neuroscientist 21, 290–305. doi: 10.1177/1073858414537560

PubMed Abstract | CrossRef Full Text | Google Scholar