Detecting local perturbations of networks in a latent hyperbolic embedding space

Any system whose behavior can be seen as resulting from the interaction of a set of interconnected units can be modeled in the mathematical form of a graph. The last two decades have seen a growing interest in applying graph theory to describe the complex structure and function of real-world interconnected systems, such as the Internet, social networks, and the brain. Evidence has revealed¹ that these systems share some peculiar characteristics, like the small-world property,² hierarchical organization,³ and efficient navigability.⁴ Therefore, understanding and analyzing them would require an expansion of standard graph models toward more specific and advanced models and methods.

We can think of the brain as a set of coarse-grained anatomical regions, connected by neuronal pathways, whose density and direction can be mapped comprehensively and non-invasively with neuroimaging techniques, namely, magnetic resonance imaging (MRI). The resulting graph encodes this information in the edge weights and represents the brain’s structure. Mapping brain functions to a graph, instead, involves the recording of brain activities, in the form of blood oxygen level (fMRI) or time-varying electrical or magnetic cortical signals (EEG, MEG). Here, nodes represent the recording sites and the edges encode some measure of correlation between signals. The resulting structural and functional brain networks are the objects of the so-called network neuroscience,⁵ whose main challenges are understanding the emergence of complex} behavior and cognition and detecting and characterizing brain diseases.⁶ 

Indeed, the presence of a neurological disease is usually reflected by alterations of local and global network properties. This has been shown for a wide variety of diseases and disorders, such as Alzheimer’s disease,^7–9 stroke,¹⁰ schizophrenia,¹¹ and epilepsy.^12,13 Comparing brain networks is a crucial task in many network neuroscience applications.¹⁴ Some diseases can be seen as a network anomaly or perturbation, and the comparison of brain data from a single patient over time (or between groups of diseased and healthy subjects) is the way one can detect and localize the onset of the disease and the starting point to predict future developments.¹⁵ In recent years, for instance, brain network analysis has been applied to improve presurgical planning and to predict the postsurgical outcome of patients.¹³ 

Network geometry, particularly hyperbolic geometry, has recently gained attention as a plausible latent geometric model of complex networks, from which typical properties, such as a heterogeneous degree distribution and strong clustering, can naturally be induced by the metric structure.¹⁶ At the same time, the approach of projecting the network onto a geometric space, called network embedding, has been used to produce highly informative low-dimensional representations of networks. Embedding networks in Euclidean spaces is a current approach to the study and comparison their connectivity structure.¹⁷ However, such embeddings often require the use of high dimensions and fail to encompass some graph topological properties, such as its hierarchical structure. In contrast, non-Euclidean mapping methods often require a lower number of dimensions to accurately embed graphs and better reflect their organization.¹⁸ 

The hypothesis is that complex network analyses based on these representations can convey additional information on the complex system at a much lower computational cost.¹⁹ When applied to structural brain networks from healthy subjects, embedding in hyperbolic space revealed a multiscale architecture of connectivity and an accurate clustering of anatomical regions in the hyperbolic disk.²⁰ However, in the case of studying neurological diseases, this method is still rarely used. A recent study on functional brain networks (obtained from magnetoencephalography data) shows its potential in diagnosing early Alzheimer’s states compared to network analysis using connectivity measures alone.²¹ However, we are not aware of any studies applying the hyperbolic embedding to epilepsy.

In this work, we study the advantage of non-Euclidean (hyperbolic) representation of networks to identify and characterize local connectivity perturbations. For this, we introduce two scores for comparing a perturbed network with its original version, the details of which are set out in Secs. II B 1 and II B 2. We assess them on the localization of single-node perturbations on synthetic networks (Sec. IV A). We finally apply them to localize disrupted areas following an epilepsy surgery (Sec. IV B).

Hyperbolic embedding is a non-linear method that takes into account the hierarchical structure of the graph (if any) and depicts both the hierarchies and the similarity of the nodes in the embedding.²² Embedding in hyperbolic space replicates properties of the graph, such as sparseness, self-similarity, scale-free degree distribution, small-world characteristics, and community structure, that embedding in Euclidean spaces cannot capture.²³ Hyperbolic embeddings have been proven to be useful for community detection and link prediction, as the embedding captures the similarity between network nodes in a graph.^23,24 Therefore, we chose hyperbolic space, instead of Euclidean one, for the network embedding. This decision allowed us to capture the properties of the original network in a lower-dimensional space and quantitatively determine local perturbations to the connectivity structure.

Hyperbolic geometry studies spaces with a constant negative curvature $K$ that do not conform to Euclidean geometry. When using hyperbolic embedding, the graph is projected onto a hyperboloid that can be further projected onto a two-dimensional hyperbolic space model like the Poincaré disk or the so-called Klein model disk. Current methods to embed networks in a hyperbolic disk essentially belong to two families:^22,25 generative model-based (e.g., Mercator²⁶) and machine learning-based. Here, we project our networks directly onto the 2D Poincaré disk using the coalescent embedding method,²⁴ which belongs to the second family, preferring it for its remarkable versatility and computational speed.

At a glance, the algorithm assigns radial coordinates to the nodes, where the radius denotes the centrality of the node, and the angle between them represents their topological proximity. In a preprocessing step, the weight of edges in the connectivity graph can be adjusted using a repulsion-attraction rule that prioritizes edges with a higher role in information transmission.²⁴ The resulting matrix is then projected onto a 2D space with a nonlinear dimension reduction method, such as Isomap or Laplacian eigenmaps.¹⁷ The nodes are adjusted equally in the polar coordinates by modifying the angular coordinates while maintaining the angular order of the nodes. The final step involves modifying each radius based on its rank, which is given by their degree. Within the disk, the angular distance between two nodes determines their level of similarity, according to a certain definition: for example, the number of neighbors they have in common. Additionally, each node’s degree centrality is proportional to its distance from the disk’s center: the closer to the center, the more central to the network. An illustration of the hyperbolic embedding of a graph is given in the left frames of Fig. 1. For the results presented in this paper, we used the dimensional reduction technique of Isomap, as it is globally optimized and preserves the global network structure by considering the shortest path distances between nodes.¹⁷ 

Proximity relations in the latent embedding space, thus, encode the importance of the interaction between nodes. In the case of the brain networks, where nodes correspond to anatomical or functional regions of interest, similar nodes that share, for instance, a high number of topological neighbors, likely belong to the same anatomical cluster and have similar characteristics as computation units in the brain. While nodes with higher centrality usually have a role of integration in the system, interacting with a larger number of anatomical regions. In the hyperbolic embedding, this is encoded in the condition that higher degree nodes are embedded to more central positions in the disk, where the hyperbolic metric naturally ensures that they are typically close, with respect to a more peripheral node, to a larger number of nodes (given that hyperbolic geodesics between peripheral points deviate toward the center).

We define the $k$-Neighborhood of node $i$ relative to the hyperbolic embedding of network $G$ as the group of $k$ nodes at a smaller hyperbolic distance from $i$⁠, and we call it $kN(i)$⁠. It represents the nodes with which node $i$ has a stronger interaction, as explained in Sec. II A.

To compute the hyperbolic distortion, we consider the same group of nodes in the embedding relative to the perturbed network $G ′$⁠, denoted as $kN( i ′)$⁠, and quantify the change in their hyperbolic distances from $i ′$ due to the perturbation. We call $i ′$ the node corresponding to $i$ in the perturbed network $G ′$⁠. We note that $kN(i)$ and $kN( i ′)$ indicate the same set of nodes, but in two different embedding configurations, meaning they will generally have different coordinates. We also specify that the set of nodes nearest to $i ′$ in the embedding space of $G ′$⁠, $k N ′( i ′)$⁠, is generally different from $kN( i ′)$⁠. These concepts are illustrated on the right side of Fig. 1.

The $kNHypD(i)$ score reflects the degree of distortion in the $k$-Neighborhood of node $i$ due to a local perturbation of its connectivity. By adding the absolute value, we model a positive contribution to the score in the case where a distance within the $k$-Neighborhood is reduced by the perturbation. This may occur due to a rewiring process: increasing the degree of a node can push it toward the center of the disk, or increasing the connectivity between two nodes can bring them closer. Dividing the difference by $dist h y p(i,j)$ gives a higher weight to nodes with which $i$ has a stronger interaction. Normalizing by a quadratic factor $dist h y p 2(i,j)$ can ensure convergence of the score for $k→N$ and $N→∞$⁠. A strong point of this measure is that it is invariant under rotations and reflections of the embedding space; therefore, no prior alignment of the two embeddings to be compared is necessary.

The best choice of $k$ in the estimation of kNHypD depends upon the data. As for other non-parametric estimators in statistics,²⁹ the variance of the estimator appears to be high for small values of $k$⁠. Although increasing the size of the neighborhood reduces the variance, it may also yield very coarse estimates of kNHypD, thus increasing the bias. Moreover, the stability gained by using higher values of $k$ can be overshadowed by the complexities involved in computing a large number of distances. Overall, based on the score’s performances across the different networks considered, it is recommended to choose a neighborhood size between 10% and 50% of the network size.

We first tested our approach on simulated networks using two models: the Watts–Strogatz model,² otherwise referred to as the small-world networks, and the Barabási–Albert model,³⁰ which consists of networks with a power-law degree distribution. The Watts–Strogatz model yields a more homogeneous connectivity, while the Barabási–Albert model produces a more heterogeneous structure and includes hubs.

We examined the change in the embedding configuration of networks under a “disconnection” perturbation, where nodes are retained (not removed, relabeled, or isolated), but some of their links are cut. This type of disruption mirrors what brain networks experience during epilepsy surgery before any connectivity rewiring mechanism begins. Specifically, we designed a random perturbation method that allows for control over the number of target nodes and the size of the perturbation. This method randomly cuts a percentage of links among those available for each target node (see Fig. 2). To ensure method compatibility with real brain networks (see below), we selected synthetic networks with similar size and density: $N=150$ nodes, parameters $k=20$ and $p=0.1$ for Watts–Strogatz, and $m=8$ for the scale-free networks. In our simulations, we excluded realizations with disconnected nodes in order to maintain a connected network.

Finally, we assessed our scores using anatomical brain connectivity networks derived from patients who underwent epilepsy surgery. We examined the brain networks of 51 patients, estimated using T1-weighted structural and diffusion-weighted MRI scans before and after surgery. The pre- and post-operative networks were obtained from an open dataset associated with the research paper by Sinha et al.¹² The dataset included 51 individuals who underwent epilepsy surgery, from which 30 patients underwent surgery for a left temporal lobe epilepsy (TLE) and 21 for a right temporal lobe epilepsy. Within each category, participants were further divided based on post-surgical results: those with positive surgical outcomes (seizure-free after surgery) and those with negative outcomes (not seizure-free after surgery). The surgical procedure involved removing brain tissue located in the presumed epileptic zone (specifically, in the temporal lobe) to prevent seizures. In graph theory, this corresponds to cutting a subset of edges within a localized section of the network. The brain networks before the surgery have the following properties: a mean degree of $28±2$ nodes, a standard deviation of the degree distribution of $12±1$⁠, and a mean clustering coefficient of $0.57±0.02$⁠. After surgery, the perturbed networks display a mean degree of $25±2$⁠, a standard deviation of the degree distribution of $11±1$⁠, and a clustering coefficient of $0.55±0.03$⁠. Brain networks have $N=114$ nodes, and there are no disconnected nodes after surgery.

When the network’s connectivity is locally disrupted, we anticipate that this will be apparent in a change in the proximity relationships of the disrupted node within its embedding space, and we expect our scores to detect this change. We hypothesize that the geometric distortion will increase statistically with the size of the disruption, both locally and globally.

To test this hypothesis, we apply the perturbation method described in Sec. III A 1 to selected synthetic connectivity models—specifically, small-world (Watts–Strogatz) and scale-free (Barabási–Albert) networks. To avoid specific perturbations from influencing the results, we repeatedly apply random perturbations to a single node a number of times that is proportional to the number of possible distinct perturbations. We estimate this number as $N p= ( k ¯ ⌊ k ¯ × m ⌉ )⋅N$⁠, where $N$ is the network’s size, and the other factor is an estimate of the possible perturbations for a node with a degree equal to the average degree $k ¯$⁠. In the binomial coefficient estimation, $m$ indicates the perturbation size as a fraction of links to cut, and $⌊x⌉$ denotes the nearest integer to $x$⁠.

We uniformly sample within the degree domain, allowing nodes with high degree $k≫ k ¯$—and, thus, a higher number of possible perturbations—to be sampled more frequently. Additionally, we choose the sample size to be a fraction $N p α$ with $α>1$⁠, ensuring both the statistical significance of the analysis and the smoothness of the computation. For a given network model (small-world or scale-free), we also average over $n G$ random network realizations.

In the trivial case where a perturbation of zero magnitude is applied, the embedding configuration of $G$ and $G ′$ is identical. As expected, both the kNHypD and HypDisp scores yield a value of $0$ uniformly across the network. Instead, the embedding configuration always changes when a non-null perturbation is applied.

To quantify how much the embedding configuration is distorted locally and globally, for each perturbation, we compute (i) the two geometric scores at the target node and (ii) the average value of the scores relative to the nearest neighbors of the target node, i.e., the nodes sharing an edge with the target node (they are just after the target node, the most physically affected by the perturbation). (iii) The baseline (BL) is defined as the average of the perturbation scores across the entire network. We then divide the results obtained for the $n p$ perturbations over $n G$ random network realizations based on the degree of the perturbed node. This enables us to observe the perturbation as a function of degree. Finally, for each degree class, we compute the average values of (i) with respect to the baseline (iii), indicated in the figure as “Target score avg−BL.” Similarly, we obtain the average value of (ii) with respect to BL, referred as “Target neighbors score avg—BL.” The measure $B L ¯$ reported in each frame indicates the value of the baseline further averaged over all the degree classes.

The results for the small-world network model are shown in Fig. 3. The HypDisp scores at the target nodes are significantly higher than the baseline scores, practically for all nodes’ degree classes. There is only a slight convergence of the target node’s score and the score from its second-order neighbors at the extremes of the degree domain, but there is never a reversal of the trend or an overshoot of the baseline score. The kNHypD score also correctly identifies a perturbation, except in the extremes of the nodes’ degree domain. If the targeted node has few connections, the perturbation cuts too few links relative to the original link density, and the score may not be sensitive enough. On the other hand, if the targeted node has more links, the perturbation results in a delocalization of the perturbation effect over many areas of the network, which can result in a target score comparable to the baseline.

The values of the average baseline $B L ¯$—defined here as the global average of the scores over all the network’s nodes and over all network realizations and perturbations—suggest that the magnitude of the perturbation definitely plays a role in the performance of the scores. As anticipated in the beginning of this section, a stronger local perturbation induces a higher distortion over larger parts of the network, which is reflected in increasing values of $B L ¯$ (this effect is reported in Sec. 2 of the supplementary material). From the bottom plots in Fig. 3, we also observe that a perturbation of higher magnitude yields a larger difference between the target kNHypD score and baseline.

The results for the scale-free network model are presented in Fig. 4. For HypDisp score, results indicate that the localization power becomes more stable as the perturbation magnitude increases (top plots in Fig. 4). Regarding the kNHypD score, there is still a notable increase in the peak distance from the baseline, but it fails to identify the perturbation when high-degree nodes or hubs are targeted. Since this model has a power-law degree distribution, perturbing a hub leads to disturbing a larger number of distant regions of the network, decreasing the difference between the error of the perturbed node and the baseline. As shown in Fig. 4, this effect appears for both scores and is more pronounced for weak perturbations.

For both network models, the scores calculated from the second-order neighbors of targeted nodes are similar to those obtained from the baseline, suggesting a good localization power, as these nodes are practically unaffected by perturbations. Strong local correlations of nodes (assortativity), however, could have an effect and induce a bias in the scores of these “second-order target nodes.” Interestingly, at the average degree, kNHypD scores show a peak in Watts–Strogatz networks and a crossing with the baseline in the Barabási–Albert networks. In contrast, no clear trend was found for the kNHypD score as a function of the target node’s clustering coefficient (results in the supplementary material). Therefore, we can conclude that (i) the proposed perturbation scores have good localization and a small bias, (ii) the performance improves with perturbation magnitude, and (iii) the performances of the kNHypD score depend on degree heterogeneity.

In the calculation of the kNHypD score, we have taken into account the displacement of the ten nearest geometric neighbors of each node, i.e., $k=10$⁠. In the supplementary material, you can find the results for two extreme values of this parameter (⁠ $k=1$ and $k=100$⁠).

To illustrate the effectiveness of our scores on real data and compare their performances, we used them on brain network connectivity perturbations due to epilepsy surgery (refer to Sec. III B for further details). For illustration purposes, here, we considered patients with a positive surgical outcome. For these patients, we compared two networks: the pre- and post-surgery networks. The network after surgery is a modified version of the original network, with fewer edges specifically removed in the surgery area. A way to assess the perturbation of node $i$ due to surgery is the perturbation degree, $P d e g(i)$⁠, defined as $P d e g(i)= ( D p r e ( i ) − D p o s t ( i ) ) D p r e ( i )$⁠, where $D p r e(i)$ and $D p o s t(i)$ denote the pre- and post-surgery degree of the $i$th node. When $P d e g(i)→0$⁠, the node’s connectivity is slightly affected by the surgery. In contrast, when $P d e g(i)→1$⁠, the node is more likely to become disconnected and significantly impacted. In both groups (left and right hemisphere surgery), 3 nodes out of 114 have a perturbation magnitude greater than $0.8$⁠, while 109 nodes have a perturbation degree below $0.4$⁠. The average perturbation magnitude across all patients and nodes is $0.118±0.032$⁠.

As our scores aim to locate the disturbed region of a network, we anticipate a correlation of our scores with the perturbation degree distribution. Results shown in Fig. 5 confirm that nodes exhibiting the highest scores (encoded by colors that tend to be clearer) correspond to the nodes exhibiting the highest perturbation degree. When surgery affects brain regions with few connections, both perturbation scores yield very low values. These findings have been consistently observed in both hemispheres. A very important point is that nodes in the surgical hemisphere clearly display a stronger perturbation score than those located on the contralateral hemisphere, independently of the perturbation degree.

The spatial mapping of brain regions’ perturbation is illustrated in Fig. 6. In both cases, one can see that nodes with the highest scores are always located in the surgically impacted region (temporal lobe). We notice, however, that the perturbation score calculated with the HypDisp score covers a large region, including all perturbed nodes, whereas the kNHypD score identifies a more reduced number of brain areas within the region directly affected by the surgery.

The characterization and localization of network perturbations have applications in various domains, such as anomaly detection (in disease networks, computer networks, etc.) and brain surgery characterization. This study introduces two scores based on hyperbolic embedding to characterize and localize local network perturbations. Results from synthetic networks validate the proficiency of the proposed scores in identifying perturbed nodes. To showcase their practical utility, the scores are applied to real brain networks from epileptic patients who underwent surgery. The results reveal that both scores successfully pinpointed the surgical area as the most perturbed region. Notably, the HypDisp score identified significant perturbations in a broader brain zone beyond the surgical region, whereas the kNHypD score highlighted perturbations in a more localized area.

The novelty of our scores resides in an explicit analysis in a latent hyperbolic space. We expect that the use of these low-dimensional and highly informative representations will be crucial for brain network studies beyond standard connectivity measures and Euclidean geometry.

For this study, the coalescent embedding method²⁴ has been selected from a range of available methods.^22,25 This method encompasses various benefits, such as (i) the absence of a minimization step that prevents errors introduced by local minima and also increases reproducibility (the same graph used as input will always produce the same embedding) and (ii) a shorter computational time with respect to methods involving a likelihood maximization, such as the Mercator.²⁶ The latter, based on a generative model, has the advantage of producing more meaningful representations. This could allow, for instance, to generate synthetic networks with similar properties to the input network. However, for the purposes of this paper, we have decided to prioritize speed and versatility over this feature. A comparison of the computational time of embeddings using coalescent and Mercator methods can be found in the supplementary material.

Similarly, we notice that other embedding techniques could potentially be combined with the proposed perturbation scores, and that alternative dimensional reduction methods (e.g., Laplacian eigenmaps¹⁷) or hyperbolic mappings (e.g., Mercator,²⁶ HyperMap,³¹ or Hydra,³² among others) might be worth considering for the representation of brain networks.³³ 

In conclusion, the article presents two scores that efficiently characterize network perturbations of the type link removal. Future work could explore the characterization of other types of perturbations, such as rewiring. Our approach may provide new insights in the study of perturbations of connected systems, including social, technological, or biological networks.

See the supplementary material for (1) the dependency of the kNHypD score on the value of parameter k, (2) the evolution of the mean baseline ( B L ¯) with the perturbation degree, (3) a comparison of computation times, as a function of network size, of the Mercator and the coalescent hyperbolic embedding method, and (4) the dependency of the kNHypD score on the clustering coefficient.

The authors are grateful to Steven Rico for useful and valuable discussions during the preparation of the study. A.L. was supported by the doctoral school EDITE from Sorbonne University. M.G. acknowledges doctoral support from the Ecole Normale Supérieure de Lyon. We thank the two anonymous referees for helpful comments.

The authors have no conflicts to disclose.

 A.L. and M.G. contributed equally to this work.

A. Longhena: Conceptualization (equal); Formal analysis (equal); Methodology (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). M. Guillemaud: Conceptualization (equal); Formal analysis (equal); Methodology (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). M. Chavez: Conceptualization (equal); Supervision (lead); Writing – review & editing (equal).

The brain network data that support the findings of this study are publicly available in the supplementary material of this article (Ref. 12).