Link Prediction based on Quantum-Inspired Ant Colony Optimization

Incomplete or partial observations of network structures pose a serious challenge to theoretical and engineering studies of real networks. To remedy the missing links in real datasets, topology-based link prediction is introduced into the studies of various networks. Due to the complexity of network structures, the accuracy and robustness of most link prediction algorithms are not satisfying enough. In this paper, we propose a quantum-inspired ant colony optimization algorithm that integrates ant colony optimization and quantum computing to predict links in networks. Extensive experiments on both synthetic and real networks show that the accuracy and robustness of the new algorithm is competitive in respect to most of the state of the art algorithms. This result suggests that the application of intelligent optimization to link prediction is promising for boosting its accuracy and robustness.

SCiEnTifiC RePORTS | (2018) 8:13389 | DOI: 10.1038/s41598-018-31254-3 robust and well-performed in general regardless of the hidden generative model of the network. This may originate from its model free property, which is adaptive since it does not assume any underlying mechanism.
In this paper, we propose a quantum-inspired ant colony optimization algorithm (QACO) to predict missing links in networks, which integrates ant colony optimization and quantum computing 28 . To sum up, our contributions are three-folds. First, based on the assumption that the pheromone of artificial ants can reflect the importance of a path, we introduce a biased ant colony optimization algorithm to remedy incomplete datasets. In our algorithm, the visibility integrates the quasi-local information of nodes. Second, quantum bits and quantum logic gates are introduced into the ant colony optimization algorithm, which can effectively prevent the optimization process from falling into a local optimum. Third, extensive experiments are conducted on many benchmark networks. Our experimental results show that the accuracy and robustness of the QACO algorithm is competitive in respect to most of the state of the art algorithms.

Results
The QACO algorithm. Overview. Consider an undirected network G = (V, E), the corresponding complete network G′ is constructed by connecting all the disconnected node pairs in G. Let all artificial ants randomly walk in G′, where links and nodes are allocated some amount of pheromone. The probability that a link or node is visited by ants is proportional to its pheromone.
Assume an ant visits n nodes, which results in a walking path. The probability that an ant travels from node v i to node v j is p ij . The value of p ij depends on the pheromone of the path, the visibility of link (v i , v j ) and the quantum pheromone of v j . After an ant reaches its destination (say v j ), the pheromone on the links and nodes in the path will be updated according to a certain rule. In turn, the updated pheromone on the links and nodes will affect the paths of the ants in the next iteration. Generally, the pheromone and visibility of the links will lead the ants to the globally optimal paths, since following the quantum pheromone is an effective way to avoid local optimums. Finally, the pheromone τ ij and visibility η ij on link (v i , v j ) can, to some extent, reflect the similarity between v i and v j .
• Quantum computation: The basic principle of quantum computation is that the relative phase and the probability amplitude of each ground state of the superposition state are not constant. The occurrence probability of each ground state varies with time, leading to the corresponding variations of the superposition state. Quantum computation is based on a device called quantum gate, which can realize a logical transformation in a certain time interval. The properties of quantum computation include interference, superposition and parallelism. In this paper, we apply quantum bits to representing the quantum pheromone of nodes and quantum rotation gate to updating quantum pheromone of ants 29,30 . • Quantum pheromone representation: In quantum computation, |0〉 and |1〉 denote two basic states of microscopic particles. An arbitrary state of a single quantum bit can be represented by a linear combination of the two basic states. The sign |〉 is called Dirac mark, which represents an eigenstate in quantum mechanics. Each quantum bit is typically in a superposition state, which is a combination of the two eigenstates. Thus, it can be represented by |ϕ〉 = α|0〉 + β|1〉, where α and β are a pair of complex numbers, which denote the probability amplitudes of quantum states. For the quantum state |ϕ〉, the probabilities of collapsing to |0〉 and |1〉 are |α| 2 and |β| 2 respectively, where |α| 2 + |β| 2 = 1. A quantum bit represents two states |0〉 and |1〉. Naturally, a quantum bit of length m can represent 2 m different states. The probability amplitude of individual j with m quantum bits is thus defined as In the QACO, the quantum pheromone is represented by the quantum bit. Let the size of population be n. The quantum pheromone is defined as p = (p 1 , p 2 , p 3 , …, p n ), where p j (j = 1, 2, …, n) is the quantum pheromone of individual j 31,32 .
• Quantum pheromone updating: The core of ant colony optimization algorithm is that the ants select their paths according to the density of pheromone, which can be applied to link prediction. For each ant, the probability of choosing a node is proportional to the density of pheromone left on the node. The update of the quantum pheromone intensity of a node can be implemented by updating the rotation angle of the quantum rotating gate in QACO. Concretely, the quantum bit of a node is updated by tuning the probability amplitude of quantum bit in the quantum rotating gate. Here, the quantum rotating gate is defined as: Here, m is the number of ants used, α β 1 is the probability amplitude of the i-th quantum bit at the t-th iteration. θ i is the rotation angle of the i-th quantum bit. Its size and orientation can be determined according to the formula θ i = Δθ * sign(α i * β i ). sign(α I * β i ) is a sign function. For α i * β i > 0, sign(α i * β i ) = 1; for α i * β i < 0, sign(α i * β i ) = −1; and for α i * β i = 0, sign(α i * β i ) = 0. Δθ usually falls in the range of [0.01π, 0.08π] 33 .

Parameters setting.
• Parameter initialization: At the beginning, we uniformly set the pheromone on all links to SCiEnTifiC RePORTS | (2018) 8:13389 | DOI:10.1038/s41598-018-31254-3 ij where δ is a constant. • On the other hand, we set the initial value of visibility of link ( where ι is a constant, varying with the topology of network. Following the settings of previous studies 34 , we set ω = 0.01. Γ(i) and Γ(j) denote the sets of i and j's neighbors, respectively. k x , k y and k z denote the degrees of x, y and z, respectively. l i→j denotes the node set in the paths from node i to node j, the length of which is 3. x and y are intermediate nodes in the set l i→j . Finally, we set the intensity of quantum pheromone at node v j to where |α j | 2 denotes the probability that the quantum state of the j-th quantum bit collapses to |0〉, while 1 − |α j | 2 denotes the probability collapsing to |1〉. • State transition rule: In each iteration, the transition probability that an ant k moves from v i to v j is defined as where τ ij and η ij are the pheromone and visibility of link (v i , v j ), respectively. μ j is the quantum pheromone intensity of node v j . The parameters λ, κ and ν are used to control the impact of the link pheromone, visibility and node pheromone intensity, respectively. • Fitness function: In each iteration, an ant walks through a path of N nodes. The path of the k-th ant is denoted by where v i ∈ V is the i-th node in the path. The fitness function is defined to evaluate each path and update the pheromone. Here, the fitness function is defined as , where C is a positive constant and d(v i ) is the degree of node v i . Generally, a path with more densely-connected nodes gets a higher score.
• Pheromone updating: After each iteration, the QACO algorithm updates the pheromone of each link according to the following formula: , ρ is the trajectory persistence, which is in the range of [0, 1). Obviously, Δτ(i, j) is proportional to the frequency that the link is visited by ants. If we consider only the fitness Q in the update of pheromone, the gaps of pheromone among links may grow. In this case, the paths of ants are likely to fall into local optimums. To remedy this drawback, β k is introduced to Eq. (6). In Eq. (6), , where the probability that ant k selects node v j is proportional to |β k | 2 . In this way, we can effectively prevent that the gap of the pheromone among links drastically grows. In other words, it can keep the optimizing process from falling into local optimums.
• Termination conditions and evaluations: Once the number of iterations exceeds a prespecified threshold N c , the optimizing process will be terminated. Finally, we can obtain a matrix τ + ε · η, which is normally referred to as the score matrix in some ref. 35 . To evaluate the performance of the QACO algorithm for link prediction, we will compare the precision of our algorithm to that of a number of existing algorithms in the performance evaluation section. The Algorithm. The proposed algorithm is outline in Algorithm 1, where the major steps are explained as follows: 1.
Step 1 (Line 1-2): In the network G of n nodes, m ants are used. For node j, its status in the t-th iteration is α The initial settings of the pheromone matrix τ and the visibility matrix η follow Eq. (2) and Eq. (3), respectively. In the experiments, we set δ, λ, κ, ν, ρ and ε to 1, 1, 2, 1, 0.9 and 0.2, respectively. ι is determined by the topology of networks, which is set to 0 or 1. 2.
Step 2 (Line 4-7): The m ants are randomly assigned to m nodes after shuffling the order of nodes. For example, to assign 3 ants to 3 of 5 nodes. We first shuffle the order of 5 nodes, say the result is 3, 1, 5, 2, 4. Then the 3 ants are assigned to the nodes of 3, 1 and 5, respectively. Step 4 (Line 20-34): After each iteration, the pheromone matrix τ of links will be updated according to Eq. (6). Similarly, the pheromone of nodes will be updated too, according to Eq. (1). 5.
Step 5 (Line 36-37): Let t = t + 1. If t < N c , then go to Step 2. Otherwise, output the matrix τ + ε · η, which will be used to calculate the precision. • Computational Complexity: Let n be the number of nodes in the network, and 〈k〉 be the average degree of nodes. Firstly, in order to calculate the visibility matrix η, all the neighbours of the second-order neighbours of each node are required to loop over. Therefore, the time complexity of calculating visibility matrix η is O(n〈k〉 3 ). Secondly, the number of qubits on each node is the same as the number of ants m. Therefore, the time complexity of initializing quantum pheromone μ is O(nm). Thirdly, the paths are the traces that the ants traverse in the network. The time complexity for forming a path in each iteration is O(n 2 ). Then the time complexity for recording m paths in N c iterations is O(N c mn 2 ). In conclusion, the overall time complexity of the QACO algorithm can be estimated as O(n〈k〉 3  Experimental results and analysis. In our experiments, we provide the performance comparison of the proposed algorithm against four state-of-the-art algorithms mentioned in Section Methods, including the algorithms based on the indices of SPM, SBM, FBM and CH. Among them, the SPM, SBM and FBM algorithms are three global algorithms for topological link prediction, while the CH algorithm is a local algorithm for topological link prediction. The recent studies have confirmed that the SPM and CH algorithms are actually one of the best-performing state-of-the-art global and local algorithms, respectively 20 . In order to verify the robustness of our algorithm, many different types of networks (the artificial, small-size, large-size and time-evolving networks) introduced in Table SI of the Supplementary Information (SI) are used. In addition, considering the high time complexity of the SBM and FBM algorithms, we only compare the SPM and CH algorithms on the large-size and time-evolving networks. Their performance will be shown as the following.
Evaluation on artificial networks. Figures 1 and 2 show the precision of the tested algorithms on the nPSO and WS networks, respectively. Figure 1 shows the average precision on the nPSO networks with 8 communities. We select the parameters (r, m, T and N) of the nPSO model, based on the topological features of the small-size real networks mentioned above. In Fig. 1, one can observe that the algorithms achieve fluctuating performances for low temperature (high clustering), relatively stable performances for medium temperature (medium clustering) and pretty stable performance for high temperature (low clustering). Furthermore, the performance of these algorithms generally decays with temperature. In general, the CH algorithm outperforms the other algorithms for most of the parameter domain. Meanwhile, as the network size N expands, the advantage of CH algorithm becomes more pronounced. Whereas, our algorithm is ranked the second. This result also confirms the robustness of our algorithm. Figure 2 shows the average precision on the Watts-Strogatz networks. The values of N and m are the same as what they are in the nPSO networks. The values of the parameter β are adjusted to build a network sharing a similar clustering coefficient and average shortest path length with the nPSO networks. In Fig. 2, one can observe that the performance of the SPM algorithm is ranked first, while that of our algorithm is ranked second again. Generally speaking, the robustness of our algorithm is better than that of the SPM and CH algorithms on both artificial networks.
Evaluation on small-size complex networks. In order to test the robustness of our algorithm, 20 small-size network datasets are collected, which covers social networks, biological networks and technology networks, etc. Table 1 shows the precision of the tested algorithms on the 20 small-size networks. One can observe that the SPM algorithm achieves the highest performance in 12 out of 20 networks, tying for the first place in two of them, whereas our algorithm receives the best result in four networks, tying for the first place in one of them. In Table 1, the mean precision of our algorithm over the tested datasets is tied with the SBM algorithm for second place, just following the SPM algorithm. The last row in Table 1 also shows the mean ranking of the tested algorithms on the networks. The precision-ranking metric is a more robust and reliable metric for assessing the overall performance. The mean ranking of the algorithms over all the networks represents the final evaluation score. One can observe that the SPM algorithm is still the best performing approach, as already deducible from mean precision, with an average ranking of 2.15. Our algorithm is ranked second with 2.78. The third approach is CH with 3.20. The SBM and FBM algorithms are ranked the fourth and fifth, respectively (The result of ranking score for each network can refer to Table SII in the SI). In short, compared with these state-of-the-art algorithms, the performance of our algorithm basically exceeds the CH, SBM and FBM algorithm, while lags behind the global algorithm, SPM. In addition, the performance of the CH algorithm (a local approach) exceeds that of the SBM and FBM algorithm (two global algorithms), which demonstrates its capability of prediction.
In order to check the statistical significance of the difference in performance between the algorithms, we perform the pairwise permutation tests (10,000 iterations), based on the precision-ranking values of each algorithm over the networks (columns of Table SII in the SI). Table SIII in the SI presents the pairwise p-values, adjusted for multiple hypothesis comparison by the Benjamini-Hochberg correction. The p-values lower than the significance level 0.05 are highlighted in bold. We find that the mean performances of the SPM-CH, SPM-SBM, SPM-FBM and QACO-FBM pairs are significantly different. In the sense, the result confirms the advantage of the SPM in performance, on the other hand, confirms that our algorithm is competitive in respect to the four state of the art algorithms. Note that the SPM is a global algorithm with time complexity, N 3 , while our algorithm is a quasi-local algorithm with time complexity, N 2 .
Comparing the performance of our algorithm in a variety of networks, we find that our algorithm is suitable for the networks with disassortative mixing 36 . For example, for the four networks, mouse neural, ACM2009 contacts, physicians innovation and haggle contacts networks, our algorithm is all ranked first in terms of precision-ranking. Their Pearson coefficients are −0.52, −0.12, −0.08 and −0.47, respectively. Conversely, when the Pearson coefficient is positive, our algorithm typically performs poorly. For example, the precision-ranking Evaluation on large-size real complex networks. Aiming to test the scalability of our algorithm, 12 large-size network datasets of different type are selected, which covers social networks, internet networks and technology networks, etc. Their topological features of common interest are likewise shown in Table SI. Table 2 shows the precision of the tested algorithms on the networks. It is evident that the CH algorithm obtains the highest performance in 6 out of 12 networks, tying the first position in two of them. The overall performance can be measured by the mean ranking. The mean ranking of the CH algorithm is 1.88, which is higher than its counterparts. With respect to the value of precision, the CH algorithm performs very poorly on the router network, which drags its mean precision. This result may be induced by the sparsity of the router network (its average clustering coefficient is 0.01), which restricts the formation of local community and ultimately leads to poor performance of CH algorithm on this network. Instead, the SPM algorithm performs very well on the yeast and arxiv astroph network. The remarkable advance leads to its mean precision is ranked first. However, its mean ranking is at the bottom, since its performance on the other networks is not desirable. The performance of the QACO algorithm is ranked second in both rankings. We also perform the pairwise permutation tests (10,000 iterations), based on the precision-ranking values of each algorithm over the large-size real networks (columns of Table SIV). From Table SV, one can find that the mean performances of the three algorithms are not significantly different on these networks. Relatively speaking, the CH algorithm performs better on most of the large-size real networks.
Evaluation on time-evolving real networks. In order to maintain the diversity of evaluation framework, we adopt the evaluation framework mentioned in previous study 20 , which considers the link-growth evolution of a real network over time. Six Autonomous systems (AS) Internet topologies collected by CAIDA are selected, which is from September 2009 to December 2010, spanning 3 months in total. Their topological features of common concern are shown in Table SI. In Table 3, one can observe that the QACO algorithm outperforms the CH and SPM algorithm, with a mean precision of 0.15. One can see that the precision grows with time, going from 0.13 to 0.16 for the QACO, from 0.11 to 0.14 for the CH and from 0.07 to 0.11 for the SPM, respectively.
In conclusion, different from the removal and re-prediction framework in which the set of missing links is artificially generated by a random procedure, here the set of links that will appear after two consecutive time points is given by ground-truth information, which makes the result even more truthful and significant, confirming the effectiveness of the QACO algorithm.
In summary, our algorithm is competitive in respect to most of the state-of-the-art algorithms in balancing precision and robustness on the tested networks.

Discussions
Inspired by ant colony optimization and quantum computing, we propose a quantum-inspired ant colony optimization algorithm for link prediction in networks. By utilizing visibility, the algorithm integrates the quasi-local structural information of individuals. By using quantum bits and quantum logic gates, it can effectively keep the optimization process from being trapped in local optimums. Compared with a series of the state-of-the-art algorithms on the artificial, small-size, large-size and time-evolving networks, our algorithm exhibits a satisfying robustness and scalability. Especially on the time-evolving real networks, our algorithm outperforms all the tested algorithms. We believe that the quantum-inspired ant colony optimization algorithm provides a new paradigm for the future studies of link prediction.
Admittedly, there may be other definitions of pheromone and visibility that can improve the performance of our algorithm. The pheromone updating strategy and the parameter selecting procedure can also be further optimized. Apart from the ant colony optimization algorithm, other intelligent optimization algorithms may be more effective on some networks. All these will be explored in our future work. Methods Link prediction problem. Consider an undirected network G(V, E) where V is a set of nodes and E is a set of links. Here, self-connections and multiple links are not considered. For each pair of disconnected nodes x, y ∈ V in the network, a link prediction algorithm assigns a score S xy , which indicates the probability of x and y connecting with each other. By sorting the scores of all disconnected node pairs, those pairs at the top of the list are more likely to be connected. In order to test the performance of the algorithm, the existing links in the network, E, are randomly divided into two sets: the training set E T and the probe set E P . Here, E = E T ∪ E P and E T ∩ E P = ∅. The algorithm estimates the scores of disconnected node pairs in G based on the information of the training set, and E P is used as the benchmark for evaluating the prediction result.
Evaluation metrics. In order to measure the accuracy of link prediction, we use the following metrics: Precision and Precision-ranking [37][38][39][40] .
• Precision refers to the fraction of correctly predicted links in the predicted links. It is defined as: where m denotes the number of correctly predicted links, L denotes the number of predicted links. For a given L, the greater precision is, the better the performance of an algorithm is. • Precision-ranking is a more robust and reliable metric for assessing the overall performance. For each network, all the algorithms are ranked by precision in descending order. The mean ranking of the algorithms over all the networks represents the final score.
Benchmark algorithms. For a comprehensive comparison, three global algorithms (SPM, SBM and FBM) and a local algorithm (CH) are considered. Let x and y be two randomly selected nodes in a network. Γ(X) and Γ(y) denote the sets of x and y's neighbors, respectively. In the following, we will briefly introduce the definitions of the algorithms mentioned above.
The SPM algorithm is a structural perturbation method that relies on a theory similar to the first-order perturbation in quantum mechanics. The idea behind this algorithm is that a missing part of the network is predictable if it does not significantly change the structural features of the observable part 26,41 . It is thus defined as: where λ k , x k and Δλ k are the eigenvalue of the observed matrix, the corresponding orthogonal normalized eigenvector and the eigenvalue of a perturbation set, respectively. The SBM algorithm is based on the assumption that the probability that two nodes are connected depends only on the groups to which they belong, and it is one of the most general network models 23 .
The FBM algorithm is a global algorithm based on the same network partitioning theory as the SBM algorithm, but it introduces a greedy strategy for an efficient sampling over the space of the possible partitions, which leads to high improvements in the computational time 25 .
The CH index is based on the assumption that two nodes are more likely to be connected if their common neighbours are densely connected 20,42 . It is thus defined as: xy CH z x y ( ) ( ) where γ(z) refers to the sub-set of z′s neighbours which are also common neighbours of x and y. |γ(z)| is defined as the local community degree of z.
Data sets. In order to validate the proposed algorithm, we provide the performance comparison over two groups of network datasets. The first group composes of two synthetic networks that are Watts-Strogatz (WS)'s small-world network 43 and Cannistraci's nonuniform popularity-similarity optimization (nPSO) network 44 . The WS small-world network is generated by rewiring r * |E| links on a regular lattice with n nodes, where r denotes the randomized rewiring probability and |E| denotes the total number of links 45,46 . The nPSO model generates synthetic networks in the hyperbolic space where heterogeneous angular node attractiveness is forced by sampling the angular coordinates from a tailored nonuniform probability distribution, and the nPSO model allows to explicitly control the size, the mixing property and the number of communities of the generated network 47 .
The second group composes of three types of real-world networks: 20 small-size, 12 large-size and 6 time-evolving real networks, which are described in detail in the SI. In addition, all real networks have been transformed into undirected, unweighted and no self-loops. Moreover, we consider only the largest component of each of the all real-world networks. The basic topological properties of the largest component in each tested network are shown in Table SI.