Communicability distance reveals hidden patterns of Alzheimer’s disease

The communicability distance between pairs of regions in human brain is used as a quantitative proxy for studying Alzheimer’s disease. Using this distance, we obtain the shortest communicability path lengths between different regions of brain networks from patients with Alzheimer’s disease (AD) and healthy cohorts (HC). We show that the shortest communicability path length is significantly better than the shortest topological path length in distinguishing AD patients from HC. Based on this approach, we identify 399 pairs of brain regions for which there are very significant changes in the shortest communicability path length after AD appears. We find that 42% of these regions interconnect both brain hemispheres, 28% connect regions inside the left hemisphere only, and 20% affect vermis connection with brain hemispheres. These findings clearly agree with the disconnection syndrome hypothesis of AD. Finally, we show that in 76.9% of damaged brain regions the shortest communicability path length drops in AD in relation to HC. This counterintuitive finding indicates that AD transforms the brain network into a more efficient system from the perspective of the transmission of the disease, because it drops the circulability of the disease factor around the brain regions in relation to its transmissibility to other regions.


INTRODUCTION
The human brain is arguably the most complex of all complex systems. At the most basic structural level of interest for neurosciences, the human brain consists of 10 11 neurons and 10 12 glial cells, which communicate through neural projections (Herculano-Houzel, 2009). These cells are then packed into local circuits (DeBello et al., 2014) or large gyri, which define anatomical and functional regions in the brain. The human brain is considered to be outstanding among mammalian brains: it is the largest than expected from body size, and it has an overdeveloped cerebral cortex representing over 80% of brain mass (González-Forero & Gardner, 2018;Herculano-Houzel, 2009). Most of the complexity of these different size scales of the human brain comes not only from the number of its components but mainly from the intricate webs of connections linking these components. The emerging field of network neuroscience studies the structural and dynamical properties of these webs observed at different size scales from a variety of noninvasive neuroimaging techniques (Bassett & Sporns, 2017). The term pathoconnectomics has been coined by Rubinov and Bullmore (2013) to describe the use of network neuroscience techniques on the analysis of abnormal brain networks (see also Stam, 2014). The goals of pathoconnectomics are not only of practical relevance as in the early diagnosis of psychiatric and developmental disorders, stroke, severe brain injury, and neurodegenerative diseases, but also in the understanding of their causal mechanisms as pointed out by Raj and Powell (2018). Due to its societal challenge, Alzheimer's disease (AD) has become a major focus of the pathoconnectomic research agenda. AD is the most common neurodegenerative disorder and it represents a major growing health problem for elderly population (Head et al., 2004;Impedovo, Pirlo, Vessio, & Angelillo, 2019;Lo et al., 2010;Rose et al., 2000). It is characterized by a continuous degradation of the patient, which starts with a preclinical stage, a phase of mild cognitive impairment (MCI), and finishes with dementia. The molecular basis of these different stages appear to be linked to the presence of β-amyloid (Aβ) in senile plaques and cerebral amyloid angiopathy, as well as tau proteins (tau) in neurofibrillary tangles (Brettschneider, Del Tredici, Lee, & Trojanowski, 2015;Jack Jr. et al., 2013;Jucker & Walker, 2011). For instance, the cognitive decline in AD correlates with tauopathy (Braak & Braak, 1991;Crystal et al., 1988;Wilcock & Esiri, 1982), while the aggregation of Aβ appears to be critical in the early stages that trigger events conducting to tauopathy, neuronal dysfunction, and dementia (Hardy & Selkoe, 2002). Then, it is plausible that these proteins, Aβ and tau, originate in a particular region of the brain and then propagate through neural fibers in a prion-like manner (Jucker & Walker, 2013;Mandelkow & Mandelkow, 1998;Walker & Jucker, 2015;Warren et al., 2013).
The hypothesis of the self-propagation of AD in combination with network neurosciences has triggered the use of epidemiological models on networks to simulate the propagation of a disease factor as AD progresses. In particular, Peraza et al. (2019)  2006; Mei, Mohagheghi, Zampieri, & Bullo, 2017;Newman, 2010), in which nodes are in two possible states, infected (I) or susceptible (S). The first corresponds to brain regions wherein the disease factor is present with high probability, while the second are those free of the disease factor but that can be infected from any infected nearest neighbor. Similar principles have guided Iturria-Medina et al. (2014) in modeling the progression of AD under the Network Diffusion Model of disease progression in dementia (Raj, Kuceyeski, & Weiner, 2012).
Here we start by adopting the SI-model for the propagation of a disease factor in AD. However, we use this model to connect with the theory of network communicability, which has Network communicability: A measure of how a perturbation in a node is transmitted to another by using all possible channels. been widely used in network neurosciences (for some applications of communicability in pathoconnectomics, see Crofts & Higham, 2009;Crofts et al., 2011;de Lange et al., 2019; Communicability distance and Alzheimer's De Reus & van den Heuvel, 2014;Johansen-Berg, Scholz, & Stagg, 2010;Li et al., 2013;Mancini et al., 2016;Qin et al., 2015;Sanchez-Rodriguez et al., 2018;Taylor et al., 2018). That is, we will provide a theoretical connection between the network communicability and the probability of a disease factor of propagating from one node to another in a network. Using this connection, we will consider a measure that accounts for the difference between the circulability of this disease factor around a given pair of nodes and its transmissibility from one region to the other. This measure is a Euclidean distance metric-communicability distance-for the corresponding pair of nodes in the network. We then find the length of the shortest communicability paths between every pair of regions in human brains for cohorts of Shortest communicability path: The shortest path between two nodes in which the edges are weighted by the communicability distance. healthy (HC) and AD individuals after appropriate normalization. We report in this work that (i) the shortest communicability path length is orders of magnitude more significant in distinguishing AD from HC than the shortest topological path length, (ii) there is a set of 399 pairs Shortest topological path: The shortest of all topological paths between two nodes.
Topological path: Sequence of nodes and edges without repetition that connect two nodes.
of regions for which there are very significant changes in the shortest communicability path length after AD, and (iii) 42% of these significant pairs of brain regions interconnect both brain hemispheres, while 28% connect regions inside the left hemisphere only, in agreement with findings related to the disconnection syndrome. Additionally, 20% of these pairs of affected Disconnection syndrome: Alteration of the connections between different brain regions rather than in the regions itself.
regions are connecting the vermis with any of the two brain hemispheres, in agreement with recent results; (iv) for 76.9% of these pairs of damaged brain regions there is an increase in the average cliquishness of the intermediate regions that connect them, which implies a significantly higher energy consumption for communication between these regions in AD than in HC.

THEORETICAL APPROACH
Here we will use indistinguishably the terms graph and network and follow the classical notation in network theory (see, for instance, Estrada, 2012c). A graph G = (V, E) is defined by a set of n nodes (vertices) V and a set of m edges E = {(u, v)|u, v ∈ V} between the nodes. The degree of a vertex is the number of edges incident to it. A walk of length k in G is a set of nodes i 1 , i 2 , . . . , i k , i k+1 such that for all 1 ≤ l ≤ k, (i l , i l+1 ) ∈ E. A closed walk is a walk for which i 1 = i k+1 . A path is a walk with no repeated nodes. The length of a path is the number of edges in that path. The shortest of all paths connecting the same pair of vertices is known as the shortest topological path. A graph is connected if there is a path connecting every pair of nodes. Here we will only consider connected graphs.
Let A be the adjacency matrix of the graph, which for simple finite graphs is symmetric, and thus its eigenvalues are real. We label the eigenvalues of A in nonincreasing order: Because the graphs considered here are connected, A is irreducible and from the Perron-Frobenius theorem we can deduce that λ 1 > λ 2 and that the leading eigenvector − → ψ 1 , which will be sometimes referred to as the Perron vector, can be chosen such that its components − → ψ 1 (u) are positive for all u ∈ V.

Susceptible-infected model
We start by considering an SI model of propagation of a disease factor as AD progresses. In this case the brain regions, represented as nodes of the graph, can be in two possible states, infected or susceptible. Susceptible brain regions are those that are free of the disease factor but that are susceptible to get infected from other regions. The infected ones are those in which the probability of disease factor is greater than zero. Let i be a node of the graph G = (V, E) and let x i (t) be the probability that node i get infected at time t from any infected nearest Communicability distance and Alzheimer's neighbor. If the infection rate is given by γ we have, (Canright & Engø-Monsen, 2006;Mei et al., 2017;Newman, 2010), It can be seen that the linearized SI model, namely ẋ(t) = γA x(t), represents an upper bound for the exact SI model, that is, which in matrix-vector is given by: with initial condition x(0) = x 0 , ∀i and ∀t. The solution of the linearized SI model ẋ ⋆ (t) = γA x ⋆ (t) is given by: However, the solution of the linearized SI model has two important problems. The first is that the semigroup e tA t≥0 is unbounded, that is, lim t→∞ e tA = ∞, which poses a major problem for the use of this linearized model as a model of the SI propagation scheme. The reason is that x(t) is a probability and as such it has to be bounded as 0 ≤ x(t) ≤ 1. The second is that the solution (4) is a good approximation to the solution of the SI model only for x ⋆ (t) ≈ 0, which makes this solution useless for following the propagation of AD.
To sort out this problem we will follow here Lee, Tenneti, and Eun (2019), who proposed the following change of variable to avoid the aforementioned problems with the solution of the linearized SI model: which is an increasing convex function. Then, as 1 − x ⋆ i is the probability that node i is not infected at a given time, the new variable y i (t) can be interpreted as the information content of the node i or surprise of not being infected (see, e.g., Cover & Thomas, 2012). Let us then suppose that at t = 0 the probability that every region of the brain gets infected is the same, that is, at the beginning every node has the same probability β to be infected and to be the one from which the disease propagates. That is, In this case, the solution of the upper bound of the SI model is: where α = 1 − β is a constant. Also, if we fix γ and t we have that this solution can be written as for constants m, ζ and b. Here the constant ζ groups the previous parameters γ, α at a given time t. We will callx i (t) = 1 − e −y i (t) the approximate solution of the SI model, which is always bounded between zero and one as needed. In order to see the differences between the exact solution of the SI model x i (t), the linearized one x ⋆ i (t) and the approximate solution after the change of variablex i (t), we plot the progression of the number of infected nodes in an Erdős-Rényi random network with n = 100 and connection probability p = 0.1. We simulate the progression of the disease starting the infection with a fraction of 0.01 infected nodes in the network. The results are illustrated in Figure 1 for infectivity rates β = 0.001 (A) and β = 0.002 (B). As can be seen, the lineralized model (dotted black line) is a bad approximation to the exact solution (broken blue line), as it quickly diverges. However, the approximate solution obtained by the change of variable (solid red line) is a tight upper bound for the exact solution, and it will be used here for further analysis.

The communicability connection
It is clear from Equation 8 that the solution y(t) of the upper bound of the SI model depends linearly of e ζ A 1. This term is the sum of the corresponding rows of the exponential of the adjacency matrix. For an individual node i it is known as the total communicability of the Exponential of the adjacency matrix: Infinite sum of powers of the adjacency matrix divided by the factorial of these powers. corresponding node (Benzi & Klymko, 2013) and it can be written as where the first term in the right-hand side is the subgraph centrality (Estrada & Rodriguez-Velazquez, 2005) of the node and the second one is the sum of the communicability functions (Estrada & Hatano, 2008) from the node i to the rest of the nodes of the network (see also Estrada, Hatano, & Benzi, 2012). In terms of the propagation of a disease factor, C i represents the circulability of the disease factor around the node i. The second term represents the transmissibility from/to the node i to/from the rest of the nodes of the network. If we concentrate on the effect of the node i on another node j, then C i represents the capacity of the node i Communicability distance and Alzheimer's of increasing the probability of infecting itself and T ij is the capacity of i of infecting j. Thus, because node j is doing the same, the term represents the difference between the capacities of both nodes of increasing the probability of infecting themself to that of infecting each other. A large value of ξ ij (ζ) indicates that the disease factor gets trapped circulating at the nodes i and j, which form two islands with little transmissibility among them. A small value, however, indicates that such transmissibility is relatively large in relation to the internal circulability at the nodes, that is, the nodes have a bridge between them. We consider that this measure is important for the study of Alzheimer's disease because it should allow us to investigate whether the disease produces a patchy environment of brain regions which form islands with little transmissibility among them. The function ξ ij (ζ) can be written as is a Euclidean distance between the nodes i and j in the network. We call it the commu-Euclidean distance: The straight-line distance between two points in a Euclidean space.

Shortest communicability paths
The communicability distance ξ ij (ζ) can be calculated for any pair of nodes (connected or not) in the graph. Thus, we can obtain a communicability distance matrix (Estrada, 2012a) where s = e ζ A 11 , e ζ A 22 , · · · , e ζ A nn is a vector whose entries are the main diagonal entries of the corresponding matrix function, 1 is an all-ones vector, and • indicates an entrywise operation. However, we assume here that in a network "information" flows through the edges of the graph, such that it uses certain paths connecting the corresponding pair of nodes. In order to find the shortest communicability paths between two nodes we proceed as follows. We construct the communicability-weighted adjacency matrix of the network : Then, the shortest communicability path between two nodes is the shortest weighted path in W. That is, the shortest communicability path between two nodes i and j for a given ζ > 0 is the path that minimizes the communicability distance between every pair of nodes in the corresponding path. We have proved analytically that when ζ → 0 the shortest communicability path between any pair of nodes i and j is identical to the shortest (topological) path between the two nodes (Silver, Akbarzadeh, & Estrada, 2018). That is, the shortest (topological) path is a special case of the shortest communicability path in a network. In this work we will consider the case ζ ≡ 1, which we will call "shortest communicability path" and the case ζ → 0 which we will call "shortest topological path." Notice that the length of the shortest communicability path is the sum of the weights (communicability distances) for the edges in that path. For an example, see Figure 2. Here we will keep ζ = 1 due to the lack of any experimental value that can guide us for selecting a more appropriate value. Also, we should have in mind that decreasing the values of this parameter close to zero will make the shortest communicability paths look very similar to shortest paths, while increasing it over unity will make these paths Figure 2. Illustration of the shortest communicability path (blue thick lines) and the shortest topological path (red thick lines) between a pair of nodes in a random geometric network. The nodes in both shortest paths are highlighted with blue color and with sizes proportional to their degrees. The rest of the nodes are in gray color and with a fixed size. The shortest topological path goes through nodes averaging node degree equal to 14. The shortest communicability path goes through nodes averaging degree 8.3. Additionally, the average subgraph centrality of the nodes in both paths are 18,837.3 (shortest topological path) and 4,636.6 (shortest communicability path), respectively. very long indeed. Thus, we left for a further work the analysis of the influence of this parameter in the study of AD.

DATASET AND IMAGE PROCESSING
The dataset used consists of diffusion-weighted imaging (DWI) scans and anatomical T1 scans Diffusion-weighted imaging (DWI): A magnetic resonance imaging (MRI) modality in which the image contrast generation is based on the diffusion of water molecules. T1 anatomical scan: Image acquired by exploiting the spin-lattice relaxation process in structural MRI. of 88 subjects, 48 healthy controls (HC) and 40 AD patients from the publicly available Alzheimer's Disease Neuroimaging Initiative (ADNI) database. After preprocessing the images, a tractography pipeline was implemented by using the MRtrix software library. The main steps of the Tractography: Image processing procedure that allows the reconstruction of fiber tracts connecting brain region pairs. whole processing, which are well established in the literature, are shown in Figure 3.
First, a denoising step was performed in order to enhance the signal-to-noise ratio of the diffusion-weighted MR signals in order to reduce the thermal noise due to the stochastic thermal motion of the water molecules and their interaction with the surrounding microstructure Communicability distance and Alzheimer's . Main steps of the image processing pipeline. DWI and T1 weighted scans are preprocessed and co-registered. Then, after the fiber orientation distribution estimation, probabilistic tractography is performed resulting in a weighted connectivity matrix. Veraart et al. (2016). Head motion and eddy current distortions were corrected by aligning the DWI images of each subject to the average b 0 image. The brain extraction tool was then used for the skull stripping of the brain (S. M. Smith, 2002). The bias-field correction was used for correcting all DWI volumes. The T1-weighted scans were preprocessed by performing the standard steps: reorientation to the standard image MNI152, automatic cropping, bias-field correction, registration to the linear and nonlinear standard space, brain extraction. The following step was the intermodal registration of the diffusion weighted and T1-weighted image.
After the preprocessing and co-registration steps, the structural connectome generation was performed. First, we generated a tissue-segmented image tailored to the anatomically constrained tractography (Zhang, Brady, & Smith, 2001). Then, we performed an unsupervised estimation of white matter, gray matter, and cerebro-spinal fluid. In the next step, the fiber orientation distribution for spherical deconvolution was estimated (Jeurissen, Tournier, Dhollander, Connelly, & Sijbers, 2014). Then a probabilistic tractography (Tournier, Calamante, & Connelly, 2010) was performed by using dynamic seeding (R. E. Smith, Tournier, Calamante, & Connelly, 2015b) and anatomically constrained tractography (R. E. Smith, Tournier, Calamante, & Connelly, 2012), which improves the tractography reconstruction by using anatomical information by means of a dynamic thresholding strategy. We applied the spherical deconvolution informed filtering of tractograms (SIFT2) methodology (R. E. Smith et al., 2015b), providing more biologically meaningful estimates of the structural connection density and a more efficient solution to the streamlines connectivity quantification problem. The obtained Streamlines: Fiber trajectories that can be reconstructed by means of tractography algorithms. streamlines were mapped through a T1 parcellation scheme by using the AAL2 atlas, (Rolls, Joliot, & Tzourio-Mazoyer, 2015), which is a revised version of the automated anatomical atlas (AAL) including 120 regions. Finally, a robust structural connectome construction was performed for generating the connectivity matrices (R. E. Smith, Tournier, Calamante, & Connelly, 2015a). The pipeline here described has been used in recent structural connectivity studies, for example, see Amico and Goñi (2018); Lella, Amoroso, Diacono, et al. (2019); Lella et al. (2020); Tipnis, Amico, Ventresca, and Goni (2018). The output was a weighted connectivity matrix for each subject. Out of these 120 nodes, 24 were removed from all networks, in order to obtain only graphs without isolated nodes. Finally, a 96 × 96 matrix M for each subject was obtained.
All matrices were binarized by considering only the edges with m ij > 0 and the adjacency matrix A was obtained. The communicability distance matrix was calculated for each binary matrix and it was multiplied by the adjacency matrix A obtaining a weighted matrix W. A shortest path algorithm was performed on this matrix, thus obtaining a matrix whose entries represent the shortest communicability paths between node pairs in the network. Starting from the adjacency matrix A, the shortest path length matrix, whose entries represent the shortest paths between node pairs, was also calculated. A group-wise statistical analysis was performed in order to find brain region pairs with a statistically significant difference between HC and AD in shortest path length and shortest communicability path length. In order to make the statistical analysis more robust, permutation tests were performed by randomly assigning subjects to the two comparison groups 1,000 times. Differences were considered significant if they did not belong to 95% of the null distribution derived from the permutation tests (corrected p value < 0.05). The false discovery rate (FDR) was used for multiple comparison correction. The same False discovery rate: The rate at which positives identified by a method are truly negative ones. study could also be done by considering the weighted matrix, but we follow here the more traditional approach on binary matrix. The weighted case could be addressed in future work.

Sensitivity analysis
Our first task here is to analyze the sensitivity of the shortest communicability and the shortest topological path lengths to detect significant changes in the brain connectivity after AD. For that purpose we proceed as follow. For each connectivity matrix we calculate both the shortest communicability path length and the shortest topological path length matrices. In this case, we use permutation tests for the statistical significance analysis. For instance, let us consider the nodes v i and v j . We then calculate the length of the communicability shortest path between these nodes for each of the healthy l H k v i , v j and diseased individuals l D k v i , v j and then obtain the respective average values,l H v i , v j andl D v i , v j . Using these values, we obtain △l v i , v j = l H v i , v j −l D v i , v j . Now we proceed to a randomization of each individual into the two classes, that is, healthy and diseased, obtaining 1,000 subsets of random HC and 1,000 subsets of random AD. We then calculate △l rand , v j are computed as before but using the random sets of healthy and diseased individuals, respectively. Finally, we compare the null distribution of △l rand v i , v j with the true value △l v i , v j . Therefore, we conclude that △l v i , v j is significant if it does not belong to 95% of the null distribution, which is carried out by calculating the p value of the permutation test. Here we consider significant the nodes with corrected p value < 0.05. We use both FDR and Bonferroni correction for the mul-Bonferroni correction: A statistical method to control false positives.
tiple comparisons correction. We do these calculations for the shortest communicability path as well as for the shortest topological path.

Effects of threshold selection and normalization
Here we first consider the effects of the thresholding process on the significance of the results obtained by using the current approach. The brain networks used in this work, as usually in many brain network studies, are constructed by using probabilistic tractography. For this reason weak connections can introduce noisy effects. Therefore, the first thing that we need to investigate is how different thresholds to transform these weighted matrices into binary (adjacency) matrices affect the results. The second important question is related to the comparison of networks with very different topological characteristics to avoid the extraction of trivial facts. That is, it is very plausible that the brain networks of AD patients differ significantly in a few "trivial" topological aspects from those of healthy individuals. For instance, the edge density can change dramatically between HC and AD networks. This may produce the false impression that AD mainly produces a sparsification of the brain network that hides important structural factors produced by the disease. To avoid these problems, we will provide a normalization of the communicability geometric parameters used in this study as described below.
First, we will proceed to change the threshold at which the adjacency matrices are generated. We start as usual by calculating the mean matrix for the HC subjects, which results in a weighted matrix in which entries range from 0 to 1. Each entry represents the frequency at which the corresponding edges occur among the HC matrices. This matrix is then thresholded by varying the threshold τ as 0 ≤ τ ≤ 0.9 obtaining a binary matrix to be used as a mask. The adjacency matrices of all subjects are then projected onto this mask. This procedure resulted in 10 sets of adjacency matrices, one set for each threshold value. We use the same thresholding procedure described in , which is also similar to the procedure used, for example, in Van Den Heuvel and Sporns (2011). In Figure 4 we illustrate the results of the statistical significance for the different thresholds studied here for both FDR (A) and Bonferroni correction (B). The first shocking result is the extremely low significance of the shortest topological paths according to both multiple comparisons correction methods for all values of the threshold. According to FDR, for almost all values of thresholds the ratio of significant node pairs is more than 35 times higher than that for the shortest topological paths, while according to the most restrictive Bonferroni correction the ratio of significant node pairs is quite dependent on the threshold value and the highest ratio of significant node pairs is obtained for τ between 0.5 and 0.9.
Then, for each value of the threshold we calculate the communicability distance matrix X (P i ) of subject i. We then proceed to normalize such matrix as follow. Let us call S (P i ) the shortest communicability path length matrix of the subject P i , let m (P i ) be the number of edges of subject i, and let A (P i ) the adjacency matrix of subject i . The average communicability distance of the edges of the network is then calculated for each subject: from which we obtain the normalized shortest communicability path length matrixŜ (P i ) as: Figure 5. Illustration of a simple graph used to explain the normalized shortest communicability path length. In red we illustrate the shortest path between the nodes 1 and 6 and in blue the shortest communicability path between the same nodes. The normalized communicability shortest path length for the path marked in blue isŜ 1,6 (G) ≈ 1.939, while for that in red isŜ 1,6 (G) ≈ 2.049.
For instance, for the network G in Figure 5 the value ofξ(G) ≈ 2.4797, then the normalized shortest communicability path length for the pair labeled as 1,6 isŜ 1,6 (G) ≈ 1.939, which corresponds to the path 1 − 9 − 8 − 7 − 6. In contrast, the normalized length of the communicability path 1 − 2 − 4 − 6 is 2.049, which clearly indicates that this is a longer path in term of the communicability distance than the path 1 − 9 − 8 − 7 − 6. Notice that 1 − 2 − 4 − 6 is the shortest topological path between the nodes 1 and 6.
We then analyze the significance for the normalized communicability distance matrices, for each threshold studied here in the two cases of the FDR and the Bonferroni. In Figure 6A  we illustrate the ratio of significant node pairs versus the threshold for FDR and in (B) the same for the Bonferroni correction. We also provide the same results for the shortest topological paths. The first interesting result is the dramatic difference between the ratios of significant node pairs obtained from the shortest communicability paths and from the topological ones. While for the shortest communicability paths we have very significant ratios for both statistical parameters for certain values of the threshold, for the shortest paths we always observe very low ratio of significant node pairs both for FDR and Bonferroni correction. The second very interesting feature of this analysis is the fact that for both statistical criteria the normalized shortest communicability paths makes a great differentiation of both groups for a threshold of τ = 0.5. Notice the nonmonotonic behavior of both statistical criteria versus the threshold, which peak at the before mentioned value. Now we focus only on the results obtained after the normalization of the communicability distance matrix and for the threshold found as providing the best results. In Figure 7 we illustrate these results using Manhattan plots for the significance of node pairs according to FDR (panels A and B) as well as for Bonferroni correction (panels C and D). Notice that the horizontal red line corresponds to the significance, that is, − ln(0.05) for FDR and − ln (0.05/κ), where κ is the number of comparisons, for Bonferroni. All the points over the red line represent significant node pairs. According to the FDR correction there are 1,551 significant node pairs for the normalized communicability distance against 120 ones according to the shortest topological paths. According to the Bonferroni correction there are 399 significant node pairs for the normalized communicability shortest paths against 71 for the topological ones. We should remark that FDR controls the expected proportion of false positives, while Bonferroni controls the overall probability of making at least one false discovery. Then, because the Bonferroni correction is a more restrictive measure, we will consider only the 399 significant node pairs identified by this measure. Notice that even with such restrictive criterion the shortest communicability path identifies more than seven times the number of significant node pairs identified by the topological shortest paths.
In closing, we observe a huge difference in the sensitivity of the communicability shortest paths with respect to the shortest topological ones to the change in the brain connectivity produced by AD. In other words, while the length of the shortest topological paths appear almost unaltered after the appearance of AD, the length of the communicability ones is affected in almost all the pairs of brain regions.
Before closing this section, we would like to remark on the huge differences produced by the normalization procedure used in the current work as a way to restrict our analysis to those highly nontrivial structural features relevant to AD. When there is no normalization in FDR case, the shortest communicability path identifies 4,524 node pairs significantly different in the two diagnostic groups, out of the 4,560 possible total node pairs. That is, 99.21% of the node pairs are affected. Conversely, the topological shortest path length only identifies 124 node pairs, which represents only 2.72% of node pairs affected. Without normalization we can also study the influence of the threshold under the ratio of significant node pairs in Bonferroni case. In this case, the best results are obtained for τ = 0.8. Then, there are 4,300 node pairs significantly different in the two diagnostic groups according to the shortest communicability path, and 95 according to the shortest topological ones. These values represent 95.61% and 2.08% node pairs affected, respectively. These results show the extraordinary value of the normalization criterion used in reducing the number of significant node pairs to a handful set of highly significant ones.

DISCOVERING STRUCTURAL PATTERNS OF ALZHEIMER'S DISEASE
Considering this threshold value, the average shortest communicability path length matrix was calculated for HC and for AD. Then, we have obtained the difference between the average matrices for the HC minus that of AD ( Figure 8A). For the sake of comparison we also obtained such differences for the shortest path lengths between every pair of nodes ( Figure 8B). We used a divergent color map centered at zero to represent these differences, and the differences for the shortest topological path length were set in the same scale of the differences of the normalized shortest communicability path length. It is interesting to note how this representation allows us to visualize a different distribution of colors in the two heat maps. In particular, if we consider the histogram of the values of one of the rows of the heat maps (the row corresponding to node 35 is considered as an example), the histogram derived from Figure 8A is a bimodal distribution ( Figure 8C), while the one derived from Figure 8B is a skewed distribution centered at zero ( Figure 8D). Moreover looking at Figure 8A, two different behaviors of the distribution of values of a single row can be observed. In particular, for some rows the average difference, between HC and AD, of normalized shortest communicability path length with the other nodes of the network is mostly positive, while for other nodes it is mostly negative. Thus, the nodes seems to be clustered according to this different behavior. In order to show how these two different groups are distributed in the brain, we have represented the brain regions on a glass brain coloring the corresponding nodes according to the median of the distribution of the values over row (Figure 9). Also the dimension of the nodes is descriptive of the median value. The cluster of nodes with the highest median values includes cerebellum, vermis, and amygdala.  From these heat maps we can observe that there are pairs of nodes for which AD increases the shortest communicability and topological paths while for others it decreases them. The difference between the distributions ofŜ(HC) andŜ(AD), that is, when P i is a healthy or AD individual, respectively, is statistically significant.
A group-wise statistical analysis using permutation tests with multiple comparison correction (both FDR and Bonferroni correction) was performed in order to find which node pairs have a significantly different value of the normalized shortest communicability path length in HC and AD. This can allow one to restrict the focus to only these node pairs, among all the possible node pairs. This procedure was applied for all thresholds.
Let us call ∆ ij =Ŝ ij (HC) −Ŝ ij (AD) the difference of the average normalized shortest path communicability distance for the edge (i, j) in both the healthy and the AD cohorts. Let us then call ∆t ij =Ŝ t ij (HC) −Ŝ t ij (AD) the difference of the average shortest topological path length. Then, we can observe that among the node pairs with statistically significant different values of the average normalized shortest communicability path length, there are both node pairs with ∆ > 0 and ∆ < 0. Instead, for the node pairs with statistically significant different values of the average shortest topological path length it is always ∆ St < 0. For example, Figure 10 shows the histograms of ∆ and ∆ St (we generally call the difference of an average shortest path measure) for the significant node pairs and for the best threshold value (τ = 0.5). Figure 10. Histogram of ∆ Sp =Ŝ p (HC) −Ŝ p (AD) for the significant node pairs for the best threshold value (τ = 0.5). Blue color refers to the difference of the average normalized shortest path communicability distance while orange color refers to difference of the average shortest topological path length. Figure 11. Glass brain visualization of the difference between the mean normalized shortest communicability path length of the significant edges in HC and AD; the edge colour is descriptive of these values. The different views show the lateral and medial sides of each hemisphere, and the dorsal and ventral side.
The values of ∆ ij for the set of 399 node pairs with significantly different average normalized shortest communicability path length for the best threshold (τ = 0.5) which were selected according to the Bonferroni correction are illustrated in Figure 11.
A deeper analysis of these node pairs provides some illuminating information about the structural influence of AD. Let us resume this information as follow. From the 399 significant node pairs considered here: 1. 110 (27.6% of significant pairs) connect regions in the left hemisphere; 2. 41 (10.3% of significant pairs) connect regions in the right hemisphere; 3. 167 (41.8 of significant pairs) connect one region of the left with one region of the right hemisphere; 4. 31 (7.8% of significant pairs) connect the vermis to the right hemisphere; 5. 50 (12.5% of significant pairs) connect the vermis to the left hemisphere.
These results indicate that almost one half of all the pairs of nodes for which there is a significant difference in the shortest communicability paths after AD connect both brain hemispheres. This result supports the disconnection hypothesis of this disease in which the damage produced by AD can be attributed not only to specific cerebral dysfunctions but also to disconnection processes between different cerebral areas (Morris, 1996;Morrison et al., 1986). In particular, the disconnection between interhemispheric regions have been widely discussed in the literature and reviewed by Delbeuck et al. (2003). They compiled evidence about the hypothesis of AD as a disconnection syndrome from neuropathological data, the electrophysiological and neuroimaging data, as well as from neuropsychological data. Some of the earliest evidence supporting this hypothesis point to the fact that there is a disorganized functional activity between the two hemispheres in the early stages of AD and a loss of positive correlations between the hemispheres, which suggest a breakdown of the interhemispheric functional association. Other evidence reviewed by Delbeuck et al. (2003) are a decrease of associative white matter fibers in the corpus callosum splenium of AD patients, the existence of a modification in the functional interactions between the hemispheres, and the existence of lower coherence between the hemispheres in mild to moderate AD patients compared to controls, suggesting a disturbance of the interhemispheric functional connectivity in AD. More recently, Wang et al. (2015) obtained experimental results that demonstrate that there are "specific patterns of interhemispheric functional connectivity changes in the AD and MCI, which can be significantly correlated with the integrity changes in the midline white matter structures." And Qiu et al. (2016) reported homotopic interhemispheric functional connectivity disruption in AD but not MCI.
Another third of these significant edges are located inside the left hemisphere. These results parallel those supporting the hypothesis that AD evolves first, faster, and more severely in the left hemisphere than in the right (Giannakopoulos, Kövari, Herrmann, Hof, & Bouras, 2009;Thompson et al., 2003). As shown by Thompson et al. (2003), the spreading waves of gray matter loss were asymmetric in both hemispheres, with the left one having significantly larger deficits and with faster local gray matter loss rates than the right one. Their finding also correlated with progressively declining cognitive status. Finally, we have found 81 pairs that connect the vermis either with the right or the left hemisphere. Recently, the role of cerebellar gray matter atrophy in AD has been studied (Mavroudis, et al., 2013;Jacobs et al., 2018). Jacobs et al. (2018) have found that in the early stages of AD the vermis and posterior lobe of the cerebellum are affected, also confirming previous results by (Mavroudis, et al., 2013), who reported severe damage in the Purkinje cells from the vermis of the cerebellum in five AD patients.

Transmissibility and circulability of AD factor
As explained before, among the pairs of nodes that have significant differences between healthy and AD cohorts, there are pairs with positive as well as with negative values of ∆ ij . Pairs of nodes for which ∆ ij > 0 corresponds to brain regions that have decreased their normalized shortest communicability path length when AD is present. Those for which ∆ ij < 0 corresponds to brain regions that have increased their normalized shortest communicability path length. A resume of our results for the pairs that increase and decrease the normalized shortest communicability path lengths is given below: 1. Left hemisphere: 31 pairs increasedŜ (P i ) and 79 pairs decreased it; 2. Right hemisphere: 6 pairs increasedŜ (P i ) and 35 pairs decreased it; 3. Left-Right hemispheres: 55 pairs increasedŜ (P i ) and 112 pairs decreased it; 4. Vermis-Left hemisphere: 0 pairs increasedŜ (P i ) and 31 pairs decreased it; 5. Vermis-Right hemisphere: 0 pairs increasedŜ (P i ) and 50 pairs decreased it.
These results indicate that in total 76.9% of all node pairs that display significant change after AD have decreased the length of the shortest communicability paths connecting them. This is a highly counterintuitive finding because it literally means that the nodes "are closer" to each other in terms of their communicability distance after the AD has appeared. In order to disentangle the meaning of this finding, we will consider the example provided in Figure 5.

Communicability distance and Alzheimer's
We will consider the removal of edges which do not disconnect the graph, which are known as cyclic edges (in the graph in Figure 5 all edges are cyclic). First we will prove the following result.
Let Γ be a graph and let Γ − e be the same graph without the cyclic edge e = {a, b}. Then, if G pq (Γ) = e γA(Γ) pq we have thatḠ pq (Γ) ≥Ḡ pq (Γ − e), whereḠ pq is the average among all pairs of nodes. The proof is given by the fact that if e = {a, b} is removed, the length of all walks between a and b will increase, and no other walk will decrease its length. Consequently, the removal of an edge in a graph will drop both the average transmissibilitȳ T (Γ) =Ḡ p =q (Γ) and the average circulabilityC (Γ) =Ḡ pp (Γ) of a disease factor. However, because the communicability distance ξ pq is the difference between the circulability around the nodes p and q, and the transmissibility between the two nodes, we have the following situations. Let ∆ξ 2 pq = ξ 2 pq (Γ) − ξ 2 pq (Γ − e), and ∆G pq = G pq (Γ) − G pq (Γ − e). Then, Consequently, when ∆ξ 2 1. ∆ξ 2 pq < 0 implies that the drop in the transmissibility of the disease factor is bigger than the drop in the circulability around the nodes. That is, in Γ − e dominates the circulability to the transmissibility compared to Γ; 2. ∆ξ 2 pq > 0 implies that the drop in the circulability of the disease factor is bigger than the drop in the transmissibility around the nodes. That is, in Γ − e dominates the transmissibility over the circulability compared to Γ.
In Table 1 we report the results that illustrate the previous reasoning for the graph in Figure 5. First, we report the change in the average shortest path length ∆L after the removal of the corresponding edges according to the node labeling in Figure 5. The removal of the edges {1, 2}, {1, 9}, and {8, 9} produce significant increase of the communicability shortest path Table 1. Values of different structural and dynamical parameters for the graph illustrated in Figure 5 to which edges have been removed. See note below. Note. The edges correspond to the labeling of the nodes in the mentioned figure. ∆L (%) is the percentage of change respect to the original graph in the average communicability shortest path length. ∆Ḡ p p (%) and ∆Ḡ pq (%) are the percentage of change with respect to the original graph for the values of G p p and G pq averaged for the nodes and edges in the shortest communicability paths. ∆t inf (%) is the time needed by a disease factor to infect all the nodes of the corresponding graph in an SI simulation by using the approximate solution described below with β = 0.005 and initial condition x i (0) = 1/9 for all i.

Communicability distance and Alzheimer's
lengths, that is, ∆ξ 2 pq < 0. The edge removals {2, 5} and {2, 4} decrease the communicability shortest paths, that is, ∆ξ 2 pq > 0. As can be seen for those graphs in which ∆ξ 2 pq < 0 the relative drop of ∆Ḡ pq , averaged for all edges in the shortest paths, is significantly bigger than that of ∆Ḡ pp , averaged for all nodes. In the case when ∆ξ 2 pq > 0, the change in ∆Ḡ pq is of the same order than ∆Ḡ pp . It is reasonable to think that edge removal increases the time needed by the disease factor to infect 100% of the nodes in the network. If we designate this time change by ∆t inf and calculate it as the minimum time at which all the nodes are infected in the graph by using the approximate solution of the SI model (upper bound), we obtain the results given in Table 1. As can be seen when ∆ξ 2 pq < 0, the increase in the infection time is dramatic, ranging from 17% to 30%. In remarkable contrast, when ∆ξ 2 pq > 0, the increase in t inf is only 3.8%. In closing, when ∆ξ 2 pq > 0 the orifinal network has been transformed into a more efficient one, from the perspective of the transmission of the disease, because in Γ − e the circulability around the nodes of the disease factor is sacrificed to the transmissibility to other nodes.
It is now important to analyze what are the consequences of the decrease in the lengths of the communicability shortest paths that we have observed for AD brains. The main implication of this observation is that AD produces some damage to the brain, not necessarily only edge removals, which somehow improves the "efficiency" of the resulting networks to propagate the disease in relation to other possible damage scenarios. If we constrain ourselves to the edge-removal scenario and use the previous graph as an example, we can conclude that AD has removed those edges which "improve" the transmissibility over the circulability of the disease factor, by affecting as least as possible the global rate of contagion in the resulting networks. This is an inference based on the analysis of data produced on AD patients and comparing them with HC. In no case does it correspond to a direct observation of this effect, and we call the attention of experimenters to try to falsify this hypothesis. We should remark that in 2013, Tomasi et al. (2013) found by using MRI "a higher degree of connectivity was associated with nonlinear increases in metabolism." Recent work in network neuroscience has proposed ways to navigate the brain avoiding the "hubs" due to their high energy consumption (Seguin, Van Den Heuvel, & Zalesky, 2018;Yan, Sporns, & Avena-Koenigsberger, 2020a). Thus, according to our findings, AD produces damage in the connectivity of the brain, which drop more significantly the cliquishness-the degree is a first order approach to it-over the transmissibility. Thus, it is also plausible that propagating the AD through the different regions of an already-damaged brain is less energetic than in the nondamaged one. All in all, the damage produced by initial stages of AD seem to improve the capacity of the disease factor to propagate through the network more so than in the undamaged one.

CONCLUSIONS
There are two main conclusions in the current work. From a theoretical perspective for network neuroscience, we have shown that the communicability function-widely used in this fieldcan be mathematically connected to the solution of an SI model of disease factor propagation in a brain network. In particular, we have shown that the communicability distance accounts for the difference between the circulability of this disease factor around a brain region (node) and its transmissibility to another region of the brain.
From an application point of view, we have provided solid evidences that the shortest communicability path length is significantly better than the shortest topological path length in distinguishing AD patients from healthy individuals. We have identified a set of 399 pairs of regions for which there are very significant changes in the shortest communicability path length after AD appears, 42% of which interconnect both brain hemispheres and 28% of which connect regions inside the left hemisphere only. These findings clearly agree with the disconnection syndrome hypothesis of AD. We have also identified 20% of affected regions that connect the vermis with any of the two brain hemispheres. The most significant result is that in 76.9% of these pairs of damaged brain regions there is an increase in the efficiency to transmit the disease factor across the brain network in relation to healthy individuals.
In closing, we hope that the current work helps to shed light on an important mechanistic aspect of AD as well as to afford a better understanding of the use of communicability functions for network neuroscience studies.