Hierarchical Aggregations for
High-Dimensional Multiplex Graph Embedding

Today’s real systems are mostly made of entities that interact with each other through multiple channels of connectivity. To better represent such complex systems, a new class of graphs, called multiplex graphs, has emerged. In a multiplex graph, nodes (entities) are connected to each other through multiple types of relations (links) [1]. We refer to these relations as dimensions. Such graphs arise in various areas, such as protein-protein interactions [2], neuroimaging [3], social networks [4], online recommendation [5]. In general, multiplexes provide a more general framework that allows rich and flexible modeling of several interconnected systems.

In recent years, many efforts have been made to design effective mining algorithms for multiplex graphs [6, 7, 8, 9, 10]. Most of these approaches are based on graph embedding techniques, which aim to project the graph dimensions to a compact and low-dimensional latent space representation. An effective embedding method should achieve exploitable representations for various prediction tasks, such as node classification [11, 12] and link prediction [13, 14]. This is not straightforward because multiplex graphs may have a large number of dimensions, each containing various complementary and / or divergent information on node interactions [7]. Therefore, an optimal embedding method should only encode the intrinsic informative structures hidden in the graph dimensions.

Over the past few years, a number of embedding methods for multiplex graphs have been proposed. Some of these methods [14, 15, 16] generalize random walks [17] to the multiplex scenario. They either perform random walks on individual dimensions and then aggregate the results or run random walks that jump from one dimension to another. However, with high-dimensional multiplex graphs, long sequences of nodes are needed to account for all dimensions. Since these methods optimize the embeddings on local patches, it becomes difficult to extract global information. Another line of work uses consensus embeddings to encode information of multiple dimensions [18, 19, 20]. A consensus embedding is the representation of a node that includes all dimensions of the graph [18]. In particular, these methods leverage regularization to force the embeddings to be similar across dimensions. This strategy is effective in some cases, but in the presence of complementary dimensions, it can hinder the quality of the final representations.

The success of deep neural networks in encoding complex data is effectively aligned with learning appropriate embeddings for multiplex graphs [21]. In this regard, Graph Neural Networks (GNNs) [22, 23] extend the deep learning framework to graph-structured data. The most successful attempts to extract latent representations for multiplexes leverage GNNs [21, 12]. Most of these methods construct separate embeddings on individual dimensions with GNNs and then aggregate the dimension-specific embeddings into a consensus one using a linear weighted summation. However, these methods consider a single and simple combination. A linear aggregation of the dimension-specific embeddings cannot capture complex relations between the graph dimensions. Moreover, it is not clear from previous work how to exploit several aggregations probably because combined representations from different sets of dimensions lead to inconsistent results [24]. In this context, we argue the existence of relevant hidden dimensions, which can be established hierarchically from the dimension-specific embeddings. Thus, a high-level relation can be seen as a composition of lower-level ones. Unfortunately, existing aggregation strategies cannot handle this compositional aspect.

Multiplex graphs are increasingly characterized by the presence of a high number of dimensions (that is, a large number of links of different types connecting multiple nodes). In this setting, informative and complex latent structures can be hidden across various dimensions. Mining such high-dimensional multiplex graphs continues to pose a challenge to existing methods. In particular, the number of possible combinations grows exponentially with the number of dimensions. Aggregating all dimensions in a single linear step can cause significant information loss. To explain these problems, we provide two illustrative examples.

Single vs Multiple: We consider the case of a co-authorship multiplex graph, where nodes are the authors of research papers, and edges indicate that two authors have co-written a paper. Edges can be scattered in multiple dimensions such that each dimension represents a publication venue (that is, two authors have co-written a paper in a specific journal or a conference). In this case, different combinations can expose different information about the authors’ research domain.

We suppose that two authors have published a paper at a data mining conference and another paper at a computer vision conference. Combining the two dimensions may suggest that the research domain of these authors revolves around “image indexing or image clustering”. If the two authors have also co-written an article in a fuzzy logic journal, then we can infer the research domain “fuzzy logic for image processing” by combining the data mining and fuzzy logic dimensions. We can see that different combinations can lead to seemingly divergent predictions. Therefore, a single combination at the initial level can cause a significant loss of information by destroying the divergent information. This information can be useful at a higher level of refinement. For example, we can infer from the initial dimension “fuzzy logic” and the hidden dimension “image indexing or image clustering” that the authors work on “robot navigation”. We illustrate this example in Figure 1(a), where nodes 1 and 4 have a high chance of sharing the aforementioned research domains. From this perspective, we can see that some initial dimensions (e.g., “fuzzy logic”) are more informative when they are used to refine higher-level dimensions (e.g.,“image clustering”).

Linear vs Non-Linear: The high-dimensionality in multiplex graphs gives rise to another phenomenon, which is the complex structure associated with information hidden in a large number of dimensions. We consider a multiplex graph that models a transportation network for a given city. Nodes represent locations in the city, and dimensions represent different means of transportation. An edge in a dimension indicates that two locations are directly linked through a transportation mean. We suppose the existence of three dimensions: Bus Lines ($G_{1}$), Metro ($G_{2}$), and Tramway ($G_{3}$). Fig. 1(b) shows an illustration of some relevant hidden dimensions that can be extracted non-linearly from the initial dimensions.

First, we can extract a relevant latent structure $G^{{}^{\prime}}_{1}$, which represents the locations that can be reached with two consecutive bus trips. This latent structure can be obtained using a non-linear operation by squaring the adjacency matrix of the Bus Lines dimension $G_{1}$. Since there is a path $2\rightarrow 1\rightarrow 4$ in $G_{1}$, the latent structure $G^{{}^{\prime}}_{1}$ contains a new link $2\rightarrow 4$. This link is not present in the initial graph structure $G_{1}$.

Another relevant latent structure $G^{{}^{\prime}}_{2}$ represents the locations that can be accessed by a metro trip followed by a tramway trip. This dimension can be formed by combining the Metro dimension $G_{2}$ and the Tramway dimension $G_{3}$ using a non-linear operation. More precisely, the links $1\rightarrow 4$ in $G_{2}$ and $4\rightarrow 2$ in $G_{3}$ form the link $1\rightarrow 2$. The link $4\rightarrow 3$ is formed by $4\rightarrow 5$ in $G_{2}$ and $5\rightarrow 3$ in $G_{3}$.

A third relevant latent structure can be revealed by combining $G^{{}^{\prime}}_{1}$ and $G^{{}^{\prime}}_{2}$. It contains locations that are linked by two bus trips, followed by a metro trip and a tramway trip. This combination forms a new path $1\rightarrow 2\rightarrow 3$. This path was not present in any of the original dimensions. Most importantly, the link $2\rightarrow 3$ cannot be obtained by any linear combination between $G_{1}$, $G_{2}$, and $G_{3}$.

Hierarchical Aggregations: High-dimensional data, such as images or audio, are well-known to be hierarchical in nature [25]. In other words, high-level features are composed of lower-level ones. For example, in image datasets, objects are combinations of motifs, which are themselves combinations of edges. Similarly to high-dimensional data, the complex structure of high-dimensional graphs inherits this compositional aspect. More precisely, there exist high-level hidden dimensions, which can be seen as non-linear compositions of some lower-level dimensions. The two real-world examples presented in Fig. 1(a) and 1(b) show the suitability of hierarchical aggregations to account for this compositional aspect.

Due to their hierarchical design, neural networks can construct compositions of low-level features and gradually generate high-level patterns [25, 26]. Although the state-of-the-art methods harness the deep learning framework to embed the initial dimensions separately, the aggregation step of these methods remains simple and can not discover the compositional structures hidden in the non-linear combinations of dimensions. To capture complex interactions between the initial dimensions, we advocate the adoption of a hierarchical aggregation method. We take inspiration from the deep learning framework to design our aggregation strategy for high-dimensional multiplex graph embedding, named HMGE (Hierarchical Multiplex Graph Embedding).

Our approach hierarchically combines the dimensions of a multiplex graph to construct new dimensions. Each layer of our model computes a more refined multiplex graph with a lower number of dimensions. The last hidden layer produces a one-dimensional graph. The node embeddings are computed by applying a standard GCN on the graph generated by the last hidden layer. Our training process leverages the hidden dimensions computed at each level to gradually extract higher-level ones. This progressive refinement allows to alleviate the information loss caused by the widely used single and linear aggregation step. Since non-linear combinations are computed from previous ones, the proposed approach can uncover complex interactions between the different relations. Finally, to encode globally relevant information present in diverse locations of the graph, we introduce mutual information maximization in the optimization module. The significance of this work can be summarized as follows.

• •

  Methodology: We propose HMGE, a novel embedding method for high-dimensional multiplex graphs. Our method relies on a hierarchical aggregation strategy to capture complex interactions between the initial dimensions. Progressive refinement of the hidden relations allows to align a large number of divergent and complementary dimensions to a consensus embedding, which in turn alleviates the information loss caused by the widely used single and linear aggregation step.

• •

  Datasets: To reflect the compositional nature of high-dimensional multiplex graphs in the experiments, we collected from various sources four multiplex graphs with an important number of dimensions. Specifically, we collected two protein-protein interaction graphs, a co-authorship graph, and a movie database graph. We made this data public for future research in the domain.

• •

  Experiments: We conducted detailed experiments on both synthetic and real high-dimensional multiplex graphs. Specifically, we evaluated HMGE on two downstream tasks: node classification and link prediction. Our results show considerable improvement compared to the state-of-the-art methods for several cases. Furthermore, our ablation study confirms that this improvement is imputed to our hierarchical aggregations.

Random walk is a popular approach to graph embedding [27], and several methods have been developed to generalize such a method to the multiplex scenario. SMNE [15] generates random walks in individual dimensions, allowing the method to learn dimension-specific node embeddings. After that, a transformation matrix is applied to obtain consensus node embeddings. MultiVERSE [14] improves SMNE by performing random walks that can jump from one dimension to another. Additionally, it supports the presence of heterogeneous nodes in the multiplex graph. PMNE [16] designs first- and second-order random walks to browse single layers. It uses the Node2Vec [28] parameters $p$ and $q$ to flexibly traverse the neighbors of the node.

A drawback of random walk-based methods is that they are inherently transductive (that is, they cannot handle the arrival of unobserved nodes). GATNE [29] solves this issue by generating training samples using random walks and then computing node embeddings with trainable transformations. In the case of high-dimensional multiplex graphs, random walk-based methods suffer from a major limitation. Indeed, to cover all the graph dimensions, long sequences of random walks must be generated. Since the skip-gram model takes advantage of local patches to optimize representations [30], it can be difficult to capture long-range dependencies in the node sequences. In other words, these methods cannot effectively uncover latent structures spanning multiple dimensions. Moreover, random walk-based methods lack the expressive power to model the compositional nature of high-dimensional multiplex graphs.

While learning node embeddings, sharing information between dimensions is essential to capture common and complementary information. To this end, a handful of methods rely on regularization. DMNE [19] uses an auto-encoder to encode each dimension of the multiplex graph. After that, a regularization term is applied to force the dimension-specific embeddings to be similar. This results in a consensus embedding for each node. MELL [20] follows a similar approach, but for directed multiplex graphs. It uses two embedding vectors for each node: one as a head vector and another as a tail vector. In addition, MELL learns layer-specific representations (i.e., embeddings for entire layers). These representations are used to compute the probability of an edge between two nodes in a specific layer. The optimized loss function is computed from these probabilities, in a similar fashion to the first-order loss of LINE [31]. MTNE [18] also uses a regularization term to encode the different dimensions to the consensus embeddings. In contrast to the other methods, the dimension-specific embeddings are generated from a pretraining stage.

Sharing information with regularization is a simple technique to build consensus embeddings. In fact, most of the regularization-based methods force the dimension-specific embeddings to be similar. However, the used similarity metrics are fixed and simple. Without prior knowledge to systematically identify the best metric, these methods generally underfit the complex interactions between the different dimensions. More precisely, these approaches cannot account for the compositional structures hidden in hierarchical combinations of the initial dimensions. Consequently, the consensus embeddings obtained based on these methods can only contain low-level and common information. In the presence of complementary dimensions, forcing the embeddings to be similar can cause significant information losses, which in turn hinders the quality of the final representations. Because of these limitations, our approach does not rely on regularization. Instead, it learns to hierarchically combine the multiplex graph dimensions to extract both common and complementary high-level information.

Learning dimension-specific node representations using GNNs and then combining these representations based on simple aggregations (e.g., linear combination, concatenation) is the gold standard for the most recent body of literature. For instance, DMGI [11] computes consensus embeddings by assigning attention weights to the node representations of each dimension. The training process of DMGI leverages InfoMax (mutual information maximization) [32] to optimize the latent representations. HDMI [12] improves on DMGI by defining a higher-order InfoMax mechanism that includes node features. SSDCM [33] applies InfoMax between local node-level representations and global cluster-aware graph summary.

A limitation of these methods is the loss of information that results from the single and linear aggregation of dimension-specific node embeddings. Different combinations of dimensions can highlight different interactions between dimensions [24]. This means that some information can only be uncovered by combining specific sets of dimensions. Second, these methods do not account for the compositional nature of high-dimensional multiplex graphs. Consequently, they cannot learn hidden structures that capture complex interactions between different dimensions. On the contrary, our approach hierarchically combines the initial dimensions, making it possible to retrieve latent structures at different levels of the hierarchy. Moreover, our method learns more complex combinations than linear aggregation methods, which we will show in Section 3.4.

Fig. 2 depicts the architecture of the proposed approach. HMGE introduces a graph neural network with $L$ layers that can tackle the multiplex case. Unlike previous methods, our approach can learn both latent features and latent graph structures. Specifically, each layer $l$ takes as input a multiplex graph $G^{(l-1)}$, whose attributes are the node embeddings $H^{(l-1)}$ ($G^{(0)}=G$ and $H^{(0)}=X$). Then, this layer outputs a new multiplex graph $G^{(l)}$, whose attributes $H^{(l)}$ and graph structures $A^{(l)}_{d}$ are computed from $G^{(l-1)}$.

Fig. 2(a) illustrates HMGE on a three-dimensional multiplex graph. The first layer computes the node embedding matrix $H^{(1)}$ and transforms the three input dimensions into two dimensions by combining the graph structures of $G_{1}$ and $G_{2}$. The second layer computes $H^{(2)}$ by refining $H^{(1)}$ using the newly created dimensions and combines the structures of $G_{1}^{(1)}$ and $G_{2}^{(1)}$ to obtain a single graph $G_{1}^{(2)}$. The final node embeddings $Z$ are computed by running a graph convolutional layer (GCN) on $G_{1}^{(2)}$.

Multiple transformations are applied inside each layer $l$, as depicted in Fig. 2(b). More precisely, each layer operates in two phases: 1) computing the node embedding matrix $H^{(l)}$; 2) generating hidden dimensions in the form of a multiplex graph $G^{(l)}$. In the next subsections, we present in detail how the embeddings are computed and how the hierarchical aggregation strategy refines the initial multiplex graph structures to obtain lower-dimensional multiplex graph structures.

The first phase of the encoding process in each layer $l$ focuses on computing the node embeddings. Given the representation $H^{(l-1)}$ and the graph structures $A^{(l-1)}_{d}$, a single graph convolutional operation is applied on each dimension $d$ of $G^{(l-1)}$ to obtain new node representations $H^{(l)}_{d}$:

The importance of each dimension in $G^{(l-1)}$ is measured by the attention weights $\beta^{(l)}_{n,d}$, which are computed as follows:

The second phase of the encoding process in each layer $l$ focuses on computing the graph structures $A^{(l)}_{d}$. As discussed in the motivation subsection, relevant information can be hidden in non-linear compositions of the graph dimensions. To identify this information, the $l^{\text{th}}$ layer computes $D_{l}$ output adjacency matrices from $D_{l-1}$ input adjacency matrices such that $D_{l-1}<D_{l}$.

Each output graph structure is computed based on a weighted summation of the input adjacency matrices followed by a non-linear activation function. The weighted summation can uncover hidden paths in the multiplex graph. It also transforms multiple adjacency matrices into a single output matrix. The activation function models non-linear interactions, and stacking multiple HMGE layers can capture more complex interactions. Accordingly, the $j^{\text{th}}$ output adjacency matrix of the $l^{\text{th}}$ layer is computed as follows:

The $l^{\text{th}}$ layer learns and computes multiple combinations, which are then used by the layer $l+1$ to improve the node embeddings $H^{(l)}$ and compute other combinations. This process forms a hierarchy of adjacency matrices as illustrated in Fig. 2(a). In addition, since each layer increments on the embeddings of the previous layer, the final representation $Z$ contains information about both the input graph $G$ and the intermediate graphs $G^{(l)}$. Moreover, by applying an activation function on the output dimensions and stacking multiple layers, HMGE can uncover non-linear latent structures hidden in the input multiplex graph dimensions.

To tackle the high-dimensionality of real-world multiplex graphs, the proposed solution gradually reduces the number of dimensions of the input at each layer. When the number of dimensions is high, aggregating all embeddings in a single step is ineffective, as it can cause a significant loss of information. HMGE addresses this issue by aggregating the embeddings multiple times on several intermediate multiplex graphs. We show in the next subsection that this hierarchical aggregation process can model more complex patterns than linear aggregation.

Most current approaches to multiplex graph embedding rely on linear aggregation to extract node representations [11], [12]. More specifically, they compute node embeddings on individual dimensions and then aggregate them using a weighted summation. Attention weights are used to account for the relevance of each dimension. In this section, we show that hierarchical aggregations can model more complex patterns than linear aggregation approaches.

For readability, we ignore the attention weights in the following equations and suppose that all graph convolutional operations have the same weight matrix $W$. Moreover, we ignore the activation functions to clearly exhibit the combinations that are formed by the different methods. For illustration purposes, consider a two-dimensional multiplex graph in the input. The following equation summarizes the embedding module of a linear aggregation method with two layers:

We can see that $A_{1}^{2}$ and $A_{2}^{2}$ are summed to compute the node embedding matrix $Z$. On the other side, a two-layer HMGE embedding module on the same input graph can be summarized as follows:

HMGE can learn more complex patterns than linear aggregation methods. Using hierarchical aggregations, both $A_{1}^{2}$ and $A_{2}^{2}$ are leveraged, while also introducing $A_{1}A_{2}$. Moreover, the weights $\alpha_{1}$ and $\alpha_{2}$ make it possible to consider the importance of each term. The term $A_{1}A_{2}$ may reveal edges that cannot be retrieved by considering $A_{1}$ and $A_{2}$ individually, or by combining them linearly. For example, as illustrated in Fig. 1(b), the combination of two adjacency matrices (specifically, in this figure, the product of the matrices) forms relevant paths between nodes that were not previously connected.

In practice, given a high-dimensional multiplex graph, several HMGE layers can be stacked to construct even more complex combinations. In the experiments section, we show that the additional expressive power of HMGE improves the performance on downstream tasks compared to linear aggregation methods on both synthetic and real-world data.

HMGE is trained based on a self-supervised strategy [35], which involves the learning of general-purpose node embeddings. Therefore, the trained model can be used to perform multiple downstream tasks. Specifically, inspired by Deep Graph Infomax (DGI) [36], we employ a loss based on mutual information maximization [32]. The training procedure maximizes the mutual information between the graph-level representation and the local patches from the set $\{z_{1},z_{2},\dots,z_{N}\}$. The graph-level representation $s$ is a summary of the entire graph $G$ and is obtained by aggregating all the node embeddings $z_{i}$:

A discriminator $\mathscr{D}$ is used as a proxy to maximize the intractable mutual information function. The discriminator takes as input a patch-summary pair and is trained to compute the probability score assigned to the pair. To this end, we select positive pairs, which consist of the graph summary $s$ and samples from $Z=\textit{HDME}(G)$. On the other hand, negative pairs are formed from $s$ and samples from $\hat{Z}=\textit{HDME}(\hat{G})$. The graph $\hat{G}$ is a corrupted version of the input graph that can be generated in multiple ways, such as swapping nodes and removing or adding links. In this work, we randomly shuffle the node features to obtain a corrupted version $\hat{X}$ and use this version as the feature matrix of $\hat{G}$. The discriminator function $\mathcal{D}$ is implemented based on bilinear scoring:

Algorithms 1, 2 and 3 summarize the proposed approach. We train the model for $T$ iterations and update the parameters $\{\alpha^{(l)},W_{d}^{(l)},V_{d}^{(l)},y_{d}^{(l)}\}$ by minimizing the loss $\mathcal{L}$ using the Adam optimizer. The time complexity of HMGE is $\mathcal{O}(TLD(\mathcal{E}M+NM^{2}+\mathcal{E}D))$, where $N$ is the number of nodes, $D$ the number of input dimensions, $L$ the number of embedding layers, $M$ the size of node embeddings, and $\mathcal{E}$ is the maximum number of edges in the graph dimensions. The memory complexity of our approach is $\mathcal{O}(LD(\mathcal{E}+NM+M^{2}+2D^{2}))$.

We first compare HMGE with DMGI and HDMI on generated data. Through this evaluation, we show that when the graph dimensions increase, the classification accuracy of the competitors rapidly decreases. In contrast, the accuracy of HMGE is significantly less affected by the increase in the number of dimensions.

In this section, we evaluate HMGE on real-world high-dimensional multiplex graphs collected from various sources (Table I). We perform two downstream tasks: link prediction and node classification.

In this section, we perform two ablation experiments to show the significance of our contributions: ablation of the whole hierarchical aggregation mechanism, and ablation of the combination weights.

Fig. 4 shows the variations of F1-macro and F1-micro scores with respect to the node embedding size $M$. We can see that the variations remain small as the scores fluctuate within the range of 71% and 82%. Since HMGE yields consistent results for a wide range of values, we can conclude that our approach is robust with respect to the node embedding size.

In this subsection, we show qualitative results to give additional insights into the inner working of HMGE.