11institutetext: School of Engineering, Stanford University, Stanford CA 94305
11email<EMAIL_ADDRESS>22institutetext: Fujitsu Laboratories of America,
Inc., Sunnyvale CA 94085
22email<EMAIL_ADDRESS>
# Ising-Based Louvain Method:
Clustering Large Graphs with Specialized Hardware

Pouya Rezazadeh Kalehbasti Work done while at Fujitsu Laboratories of America,
Inc11 Hayato Ushijima-Mwesigwa 22 Avradip Mandal 22 Indradeep Ghosh 22

###### Abstract

Recent advances in specialized hardware for solving optimization problems such
quantum computers, quantum annealers, and CMOS annealers give rise to new ways
for solving real-word complex problems. However, given current and near-term
hardware limitations, the number of variables required to express a large
real-world problem easily exceeds the hardware capabilities, thus hybrid
methods are usually developed in order to utilize the hardware. In this work,
we advocate for the development of hybrid methods that are built on top of the
frameworks of existing state-of-art heuristics, thereby improving these
methods. We demonstrate this by building on the so called Louvain method,
which is one of the most popular algorithms for the Community detection
problem and develop and Ising-based Louvain method. The proposed method
outperforms two state-of-the-art community detection algorithms in clustering
several small to large-scale graphs. The results show promise in adapting the
same optimization approach to other unsupervised learning heuristics to
improve their performance.

###### Keywords:

Graphs Community detection Clustering Ising Model

##  1 Introduction

Graphs are powerful tools for modeling and analyzing complex systems from
social networks to biological structures. In graphs, some groups of nodes can
show modular behavior [16] in that nodes in these groups are highly
interrelated and show similar traits and strong connections [6]. Communities
or clusters are a form of these modules in graphs which are identified as
clusters of nodes with a higher density of internal edges compared to external
edges [24, 8]. These highly inter-connected nodes reveal underlying structures
within a network ranging from friend groups in social networks [22], similar
proteins in biochemical structures [14], and papers of similar categories in
citation graphs [30]. _Community Detection_ (CD) is an unsupervised learning
method which aims to identify these interwoven structures (communities, aka
clusters) in graphs. CD is an NP-Hard combinatorial optimization problem.

In recent years, we have seen the emergence of specialized hardware designed
for solving combinatorial optimization problem. Examples of these novel
computing hardware include, adiabatic quantum computers, CMOS annealers,
memristive circuits, and optical parametric oscillators, that are designed to
solve optimization problems formulated as an Ising or QUBO mathematical model.
Given that many well-known problems in graphs can easily be model in this
form, there has been a growing interest in formulating and evaluating these
problem and their subsequent QUBO models on the different specialized hardware
platforms [15, 32, 34, 4, 21]. However, many real-world problems require
significantly more variables than these devices can handle, thus hybrid
methods are usually used.

##### Our contribution:

This paper introduces a hybrid method, inspired by the Louvain algorithm, for
solving the community detection problem on specialized hardware. Following the
work on leveraging specialized hardware for local search [18], the proposed
method he proposed method is expected to produce more optimal clusterings
(i.e. communities) compared to the original algorithm. This improvement is
expected since instead of using a greedy approach for assigning nodes to
communities iteratively, this method creates local optimization problems for a
set of multiple nodes and multiple candidate clusters concurrently. The method
is used to cluster several benchmark graphs where the results will show
improved modularity scores compared to Louvain algorithm across several runs.

###  1.1 Related Work

Researchers have developed several methods for finding high-quality
communities that can reveal the interworking and internal structures of
complex networks. Increasing availability of large-scale graphs, like the
World Wide Web, or social networks of Twitter, Facebook, and Instagram with
millions or billions of nodes and edges, have made traditional approaches for
community detection unwieldy in terms of computational and spatial complexity
[3]. Consequently, researchers have devised efficient algorithms, $e.g.$
Girvan-Newman [8], Louvain [3], and Leiden [33] algorithms, for effectively
detecting communities within large graphs. Some of the other well-known CD
algorithms listed by Javed et al. [11] include Hierarchical Clustering [12],
Graph Partitioning [23], Spectral Optimization [26], Simulated Annealing [9],
Clique Percolation Method (CPM) [27], and Fuzzy Detection [29]. Many of these
algorithm use _modularity_ as a metric for solving the optimization problem of
CD: such algorithms try to assign nodes to clusters which can maximize the
modularity of the entire graph [26, 23].

##  2 Methods

###  2.1 Modularity Maximization

The objective function selected here for clustering is Modularity which was
introduced by Newman [26]. This metric favors highly interconnected nodes to
be clustered together. For a graph, $G=(V,E)$, with $|V|$ the set of nodes and
edge-set $E$, the Modularity Maximization problem maximizes the following
objective function:

$Q=\frac{1}{2m}\sum_{i,j}(A_{i,j}-\frac{k_{i}k_{j}}{2m})\delta(c_{i},c_{j})$
(1)

here, $A_{i,j}$ is the $(i,j)$ element in the adjacency matrix of $G$. $k_{i}$
is the weighted degree of node $i\in V$, and $m=|E|$. In addition, $c_{i}$ is
the community node $i$ belongs to, and $\delta(u,v)$ equals 1 if $u=v$, and
equals 0 otherwise.

###  2.2 QUBO Models

Quadratic Unconstrained Binary Optimization (QUBO) models include 0/1 binary
variables $q_{i}$, biases $c_{i}$, and couplers $c_{i,j}$. The objective is to
minimize the following:

$E(q_{1},\dots,q_{n})=\sum_{i=1}^{n}c_{i}q_{i}+\sum_{i<j}c_{i,j}q_{i}q_{j}$
(2)

QUBO models can be converted to Ising models (variables are $\pm 1$) and vice
versa [2]. The objective function of the Modularity Maximization for at most
2-clusters is naturally defined in the Ising Model [21]. In addition, the
problem of at most $k$ clusters forms a quadratic function, thus leads to a
natural formulation as a QUBO. Subsequently, there has been a large body of
research utilizing specialized hardware for tackling this problem ranging from
both small to large scale graphs [34, 31, 32, 21, 35].

###  2.3 Solving Community Detection on the Fujitsu Digital Annealer

The Fujitsu Digital Annealer (DA) is a specialized architecture for
combinatorial optimization problems formulated as QUBO [1, 7]. We use the
second generation of the DA that is capable of solving problems with up to
8192 variables and up to 64 bits of precision. The DA has previously been used
in different areas such as communication [20], signal processing [28], and
data mining [4, 5].

###  2.4 Baseline Method: Louvain Algorithm

Louvain Algorithm (aka Louvain) is a hierarchical clustering method that uses,
mainly, modularity for clustering large graphs [3, 25]. Louvain identifies the
optimal clusters for nodes in an agglomerative hierarchical approach: it first
assigns each node to an individual cluster, then identifies the optimal
cluster for each node in its immediate neighboring clusters, and after several
passes over the graph, it aggregates the nodes belonging to each cluster into
super-nodes [3]. The same process is applied to these super-nodes until all
clusters become singletons and no nodes can be assigned to clusters other than
their own. The resulting super-nodes are the ultimate clusters for their
associated sub-nodes. Blondel et al. [3] and Traag et al. [33] provide
detailed illustrations of the Louvain algorithm. Next section introduces the
Ising-Louvain where an optimization-based approach replaces the greedy
approach in the original Louvain.

###  2.5 Ising-Louvain Method

The proposed method in this paper, Ising-Louvain, is an agglomerative
hierarchical clustering method. This method can be summarized in the following
steps

  1. 1.

Initially, each node is assigned to a separate cluster

  2. 2.

A breadth first Search is performed for each node to find its neighboring
nodes, and selects a set of candidate clusters for each of these selected
nodes

  3. 3.

These candidate clusters get assigned to the selected nodes in a way that the
new assignments maximize the modularity of the entire graph

  4. 4.

Steps 2 and 3 are repeated until no more improvement in the modularity of the
graph is possible

  5. 5.

Afterwards, the nodes belonging to the same cluster get aggregated into super-
nodes, and steps 2-4 are applied to these super-nodes

  6. 6.

Steps 2-5 get repeated until all clusters become singletons.

In the steps above, Step 3 is where the QUBOs are formed and solved.
Specialized hardware, e.g. Fujitsu’s Digital Annealer, can effectively solve
such optimization problems and yield high-quality solutions.

Figure 1 shows a case where Louvain cannot improve the clustering of a graph
in its current situation, while Ising-Louvain can. Here, Ising-Louvain can
find a more optimal clustering by considering two nodes [or more than one node
in general] at a time and changing the clusters of both nodes from violet to
yellow. However, Louvain only considers a single node at each step and finds
no gain in modularity in assigning a new cluster to that single node, and
hence can never achieve the more optimal clustering that Ising-Louvain can
obtain.

Figure 1: Sample case where only Ising-Louvain can improve the graph while
Louvain stalls

###  2.6 Formulating Ising-Louvain

####  2.6.1 Deriving the Objective Function

What follows formulates the modularity as the optimization objective for the
QUBO used in Ising-Louvain algorithm. Assume we have a graph with $n$ nodes
and an edge set $E$. We start by selecting a set $S$ of nodes as our ‘free
nodes’ and by checking if they should belong to a specific community $C$ or
not. The community-membership of the rest of the nodes are assumed to be
constant. Let $\boldsymbol{x}\in\\{0,1\\}$ be a vector of $n$ binary variables
indicating whether or not each node in the graph belongs to community $C$,
such that $x_{i}=0$ if $i\in C$ and $x_{i}=0$ otherwise. Following a similar
derivation as Negre et al. [21], the objective function defined in Equation 1
can be written as follows:

$\max(Q)\equiv\min_{\boldsymbol{x}}\Big{(}\sum_{i\in S}\sum_{j\in
S}x_{i}(\frac{k_{i}k_{j}}{2m}-A_{i,j})x_{j}+\sum_{i\in
S}x_{i}\big{(}2\sum_{j\in C}(\frac{k_{i}k_{j}}{2m}-A_{i,j})\big{)}\Big{)}$ (3)

We can apply the same process to nodes in set $S_{l}$ and any other cluster
$C_{l}$ and obtain similar phrases to Equation 3. When considering more than
one cluster in our optimization problem, we should consider the constraint
that each node in $S$ (defined as $\bigcup_{l}S_{l}$) should belong to one and
only one cluster at each instant. Equation 4 shows the objective function
resulting from superposing the objective functions for a set of candidate
clusters, called set $L$, with a penalty term imposing the mentioned
constraint:

$\begin{aligned} \min_{\boldsymbol{x}}\bigg{(}&\sum_{l\in L}\Big{(}\sum_{i\in
S_{l}}\sum_{j\in
S_{l}}x_{i,l}(\frac{k_{i}k_{j}}{2m}-A_{i,j})x_{j,l}+\sum_{i\in
S_{l}}x_{i,l}\big{(}2\sum_{j\in
C_{l}}(\frac{k_{i}k_{j}}{2m}-A_{i,j})\big{)}\Big{)}\\\ &+\gamma\sum_{i\in
S_{l}}(\sum_{l\in L}x_{i,l}-1)^{2}\bigg{)}\end{aligned}$ (4)

where $l$ indicates the cluster number, $S_{l}$ indicates the nodes for which
cluster $C_{l}$ is a candidate, $x_{i,l}$ is the i’th element of vector
$x_{l}$ which contains the binary variables indicating the membership of all
nodes to cluster $C_{l}$. The last term in this equation imposes the
constraint mentioned above with $\gamma$ being the penalty coefficient.

Equation 4 has a space complexity of
$O\big{(}\max{\\{|S|^{2},\mathscr{L}\\}}\big{)}$ where
$\mathscr{L}=\sum_{l}|S_{l}|$. The optimization objective has $\mathscr{L}$
variables, so the QUBO is of size $\mathscr{L}\times\mathscr{L}$.

####  2.6.2 Pseudocode

Algorithm 1 shows the pseudocode of the proposed method [partially following
the pseudocodes of Blondel et al. [3] and Traag et al. [33]]. Table 1
describes functions and auxiliary variables used in Algorithm 1, and Table 2
describes the hyper-parameters defined to fine-tune the behavior of the
algorithm.

Table 1: Definition of functions and auxilliary variables used in Algorithm 1 Function/Variable Name | Definition
---|---
Singleton(G) | Create singleton clusters for graph G
Aggregate(G, P) | Aggregate nodes in the same clusters (based on partitioning P) in
| graph G into supernodes
modularity(G, P) | Modularity of graph G with partitioning P
ModGain(P, X) | Modularity gain if partition P is updated with assignments X
modified | flag indicating if graph is modified throughout the inner loop
done | flag indicating that graph can’t be improved anymore, i.e. all clusters
| have stayed as singletons after running the inner while-loop
$nodes$ | list of nodes in graph G
Table 2: Definition of hyper-parameters used in Algorithm 1 Hyperparameter Name | Definition
---|---
Max_Nodes | Maximum number of nodes that the algorithm selects as set $S$
Max_Clusters | Maximum number of candidate clusters that the algorithm considers
| for each node in set $S$
Max_Node_Visits | Maximum number of times the algorithm includes each node of the
| graph in set $S$ over a single pass of the graph during the inner while-loop
random_seed | Seed number used to initialize the random generator in Numpy
solver_timeout | Timeout limit for each call to the solver
BFS_Depth | Depth of the BFS run by the algorithm (explained more in Section 2.7.1)
gamma | the penalty coefficient in QUBO
$counter\\_max_{out}$ | maximum allowed number of iterations for the outer while-loop
$counter\\_max_{in}$ | maximum allowed number of iterations for the inner while-loop
$\theta$ | minimum acceptable improvement in modularity over an inner or outer
| while-loop

Function  _IsingLouvain(_Graph G_)_:

P = Singleton(G);

P = RefineAndCoarsen(G, P);

return P

Function  _RefineAndCoarsen(_Graph G, Partition P_)_:

$counter_{out}=0$;

$\Delta Q_{out}=\theta$;

$done=False$;

/* Outer While Loop */

while  _(($done==False$) & ($counter_{out}\leq counter\\_max_{out}$) &
($\Delta Q_{out}\geq\theta$))_ do

$counter_{out}\mathrel{+}=1$;

$counter_{in}=0$;

$Q_{old,out}=modularity(G,P)$;

$\Delta Q_{in}=\theta$;

$modified=True$;

/* Inner While Loop */

while  _(($modified==True$) & ($counter\leq counter\\_max$) & ($\Delta
Q_{in}\geq\theta$))_ do

$counter_{in}\mathrel{+}=1$;

$modified=False$;

$Q_{old,in}=modularity(G,P)$;

/* Doing a Single Pass */

P, modified = RunOnePass(G, P, modified);

$Q_{new,in}=modularity(G,P)$;

$\Delta Q_{in}=Q_{new,in}-Q_{old,in}$;

$Q_{old,in}=Q_{new,in}$;

end while

$Q_{new,out}=modularity(G,P)$;

$\Delta Q_{out}=Q_{new,out}-Q_{old,out}$;

$Q_{old,out}=Q_{new,out}$;

if  _$|P|==n$_ then

done $\leftarrow$ True;

else

G = Aggregate(G, P);

P = Singleton(G);

end if

end while

return G, P

Function  _RunOnePass(_Graph G, Partition P, Bool mod_)_:

for  _$i\in nodes$_ do

Select $S$ and $L$;

Calculate
$\boldsymbol{\mathscr{B}}_{\mathscr{L}\times\mathscr{L}},\boldsymbol{\mathscr{B}}_{\mathscr{C}\times\mathscr{L}},\boldsymbol{\mathscr{Z}},\boldsymbol{\zeta}^{T}$;

$\boldsymbol{X}=[\boldsymbol{x}_{l,S_{l}}]_{\mathscr{L}\times 1}$ (i.e.
solution to Equation 4);

if  _$ModGain(P,X) >0$_ then

Update P with $X$;

mod $\leftarrow$ True;

end if

end for

return P, mod

Algorithm 1 Ising-Louvain

###  2.7 Implementation Details of Ising-Louvain

The first part of this section will study the strategies tried for selecting
candidate nodes and clusters and the ones which yielded the best results, the
second part studies the rule-of-thumb devised for tuning the $\gamma$
hyperparameter in Equation (4). Lastly, we summarise the different approaches
taken to reduce the time complexity of the proposed method.

####  2.7.1 Node and Cluster Selection Strategies

The authors noticed that the strategy for selecting nodes and candidate
clusters has a sizable impact on the run-time and performance (measured by the
obtained optimal modularity value) of the devised clustering algorithm. The
explored node-selection-strategies for forming set $S$ were as follows:

  1. 1.

Random strategy: selecting set $S$ at random from the graph

  2. 2.

Sliding-window strategy: selecting set $S$ from consecutive nodes in the
sorted list of nodes, and shifting the list by a random number after each pass
over the entire graph [shifting would result in (probably) different
selections during each pass]

  3. 3.

BFS strategy: selecting set $S$ based on a breadth-first-search (BFS), of
specified depth, originating from each node selected sequentially from the
shuffled list of nodes

Experiments showed that the third strategy yielded consistently better results
across different graphs and across different runs on the same graphs. The
experiments reported later on are conducted with this approach for node
selection.

The following strategies, for selecting the cluster-set $L$, were considered
in this work.

  1. 1.

BFS strategy: selecting set $L$ at random from clusters resulting from a BFS
[of specified depth] originating in each node, while assuring the current
clusters for nodes in set $S$ are included in the candidate clusters for each
node [this extra step enabled obtaining a trivial solution to the QUBO with
all nodes kept in their current clusters]

  2. 2.

Semi-greedy strategy: composing set $L$ from the top $k$ clusters for each
node in set $S$ based on a greedy search for that node [the greedy search
would inspect all neighboring clusters as well as the original cluster for the
node, and would select the top $k$ clusters which resulted in the highest
gains in modularity if the node was to be moved to those clusters]

The second strategy resulted in consistently better results across different
graphs and across different runs on the same graphs, so this strategy was used
to conduct the experiments reported later in the paper.

####  2.7.2 Strategy for Tuning $\gamma$

The authors devised the following rule-of-thumb for setting the initial value
of $\gamma$ in Equation 4:

$\gamma=max(k_{i})\ \text{for}\ i\in G$

This initial value of $\gamma$ resulted in consistent satisfaction of the
constraint imposed on the QUBO in Equation 4. Also, after aggregating the
graph at the end of each inner while-loop (Algorithm 1), $\gamma$ is updated
to account for the increasing degrees of nodes:

$\gamma_{new}=max(k_{i})\ \text{for}\ i\in G_{aggregated}$ (5)

####  2.7.3 Strategies for Reducing Runtime

The majority of the different types of specialized hardware are mainly
accessible via the cloud, commonly provided by some WebAPI access to the
hardware. In order to make proposed approach applicable to a variety of
specialized-hardware platforms, we take this into account and outline
strategies that minimize the number of calls to the hardware, and if possible,
minimize the size of the QUBO per iteration.

  1. 1.

When set $S$ had only a single node, a greedy search replaced the call to the
solver, and the node was assigned to the top choice determined by the greedy
search. This strategies is especially helpful in reducing the number of calls
to the solver when many nodes in the graph have been visited enough number of
times by the algorithm, and hence the algorithm frequently selects single-
member $S$ sets.

  2. 2.

When creating the candidate clusters set ($L$), nodes with only a single
candidate cluster were excluded from set $S$, and instead were assigned to
their single candidate clusters and were treated as normal nodes with
determined clusters. After solving the QUBO for the reduced set $S$, the
assignments for the single-candidate nodes were added back to the assignments
determined by the solver. This strategy led to up to 50-60% reductions in the
size of the QUBOs sent to the solver for large QUBOs, and up to 30-40%
reduction in the number of calls made to the solver. This strategy yields more
highlighted results after the first few passes over the graph in the inner
while-loop (see Algorithm 1), when the number of clusters in the graph shrinks
significantly, and consequently the algorithm identifies many single-candidate
nodes.

##  3 Results and Discussion

###  3.1 Experiments

Table 4 lists the 12 benchmark graphs that are used in this paper: [Zachary’s]
Karate Club [37], Meredith [19], Les Miserables [13], and the rest of the
graphs were selected from SNAP library [17]. The original Louvain algorithm
was used to produce the baseline results of clustering, and Leiden algorithm
[33] (an updated version of Louvain) was used to produce rivaling results to
compare the performance of the Ising-Louvain method against. To obtain the
Louvain and Leiden results, each were run 20-30 times and the best results
were selected. To obtain the Ising-Louvain results, several experiments were
conducted to find the best combination of hyperparameters, and the results of
the best experiments were reported. The Ising-Louvain algorithm was developed
in Python 3.7 [36] using NetworkX v2.4 [10] as the main library for handling
the graphs, and the QUBOs were solved using Fujitsu’s Digital Annealer.

In Table 4, #Clusters columns show the final number of clusters obtained by
each method, Avg QUBO Size shows the average number of variables in QUBOs sent
to the solver, and #Solver-Calls shows the number of calls to the solver.

###  3.2 Results

Table 4 shows that in all cases, Ising-Louvain has been able to obtain a more
optimal clustering than that obtained by Louvain. The only exceptions are the
three smallest graphs (Karate Club, Meredith, and Les Miserables) on which
both Louvain and Ising-Louvain have obtained the globally optimal clustering;
thus, the final modularities obtained by both methods are equal. This is
largely because Ising-Louvain, using local optimization, can consider and
evaluate orders of magnitude more possible assignments of clusters to nodes
than Louvain is able to consider with a greedy approach. This way, Ising-
Louvain can find more optimal clusterings than the original Louvain.

Figure 2: Boxplot of sizes of QUBOs solved for clustering the Facebook graph

For a more detailed analysis of QUBO sizes, Figure 2 shows the box-plot of the
sizes of all QUBOs solved when clustering the Facebook graph (c.f. Experiment
4 in Table 4): the box starts from the first quantile (Q1 = 4) and continues
until the third quantile (Q3 = 14) of the sizes with a blue line indicating
the median value (Q2 = 6). The lower and upper whiskers correspond to,
respectively, $Q1-1.5*(Q3-Q1)$ and $Q3+1.5*(Q3-Q1)$, and the circles show the
data lying outside the boundary of the whiskers. Figure 2 shows that 75% of
the QUBOs sent to the solver had less than 14 variables, and most of the QUBOs
had less than 50 variables.

Table 3: Benchmark instances Graph Name | KarateClub | Meredith | LesMiserables | Facebook | Autonomous | LastFM
---|---|---|---|---|---|---
$|V|$ | 34 | 70 | 77 | 4,039 | 6,474 | 7,624
$|E|$ | 78 | 140 | 254 | 88,234 | 13,895 | 27,806
Graph Name | arXiv | arXiv2 | AstroPh | Enron | DBLP | Amazon
$|V|$ | 9,877 | 12,008 | 18,772 | 36,692 | 317,080 | 334,863
$|E|$ | 25,998 | 118,521 | 198,110 | 183,831 | 1,049,866 | 925,872
Table 4: Summary of the results of running Louvain, Leiden, and Ising-Louvain on the benchmark graphs, along with statistics of the solver | Optimal Modularity Results | Solver Stats
---|---|---
Graph | Leiden | Louvain | Ising-Louvain |  Avg QUBO Size | #Solver-Calls
KarateClub | 0.4198 | 0.4198 | 0.4198 | 75 | 4
Meredith | 0.7457 | 0.7571 | 0.7571 | 164 | 4
LesMiserables | 0.5600 | 0.5667 | 0.5667 | 85 | 14
Facebook | 0.8356 | 0.8350 | 0.8358 | 19 | 584
Autonomous | 0.6662 | 0.6553 | 0.6572 | 7 | 1626
LastFM | 0.8170 | 0.8155 | 0.8168 | 16 | 19036
arXiv | 0.7759 | 0.7699 | 0.7756 | 11 | 36373
arXiv2 | 0.6649 | 0.6622 | 0.6650 | 32 | 115898
AstroPh | 0.6367 | 0.6309 | 0.6378 | 44 | 351410
Enron | 0.6265 | 0.6204 | 0.6285 | 23 | 146451
DBLP | 0.8301 | 0.8220 | 0.8243 | 13 | 8369823
Amazon | 0.9309 | 0.9263 | 0.9277 | 14 | 7812919

##  4 Conclusion

This paper introduces a new clustering algorithm based on the popular Louvain
clustering algorithm for large-scale graphs using quadratic unconstrained
binary optimization (QUBO). This method, called Ising-Louvain, is particularly
fit for being used with specialized hardware, including Digital Annealers,
which are capable of solving combinatorial optimization problems. The method
builds around replacing with local optimization a greedy approach taken by
most state-of-the-art clustering algorithms. The method proves effective in
obtaining improved results (as measured by modularity metric) compared to
Louvain and Leiden methods — the state-of-the-art in large-scale graph
clustering — when tested on 12 graphs with sizes ranging from 34 to 335k
nodes. The improvement is attributed to the ability of the proposed method in
evaluating orders-of-magnitude more possibilities of assigning clusters to
individual nodes, and hence being able to find more optimal clusterings than
the greedy method can find.

One noteworthy extension to this work can focus on using greedy approaches,
like Leiden and Louvain, to run the initial iterations of the proposed
algorithm on very large graphs, when the number of neighboring clusters for
each node are large and hence the created QUBOs are massive. This can save a
lot of computation time and can allow more experiments to be run on these
graphs for obtaining more optimal solutions.

Future work can also try fitting the proposed algorithm with optimization
objectives other than modularity, such as betweenness, to see the behavior of
the algorithm with those objective functions. Further, the distance metric
used in this paper was the adjacency matrix. Future research can try focusing
on graphs in which other measures of distance between nodes, such as Euclidean
distance, can be applied for selecting candidate nodes and clusters and
forming the QUBOs. This can extend the application of Ising-Louvain to,
especially, Geospatial Clustering.

Finally, the success of this paper in revising and improving Louvain algorithm
shows promise for applying the same local-optimization-oriented approach to
other unsupervised learning algorithms which depend on local search
heuristics.

## References

  * [1] Aramon, M., Rosenberg, G., Valiante, E., Miyazawa, T., Tamura, H., Katzgraber, H.G.: Physics-inspired optimization for quadratic unconstrained problems using a digital annealer. Frontiers in Physics 7, 48 (2019)
  * [2] Bian, Z., Chudak, F., Macready, W.G., Rose, G.: The ising model: teaching an old problem new tricks. D-wave systems 2 (2010)
  * [3] Blondel, V.D., Guillaume, J.L., Lambiotte, R., Lefebvre, E.: Fast unfolding of communities in large networks. Journal of statistical mechanics: theory and experiment 2008(10), P10008 (2008)
  * [4] Cohen, E., Mandal, A., Ushijima-Mwesigwa, H., Roy, A.: Ising-based consensus clustering on specialized hardware. In: IDA. pp. 106–118. Springer (2020)
  * [5] Cohen, E., Senderovich, A., Beck, J.C.: An ising framework for constrained clustering on special purpose hardware. In: CPAIOR. pp. 130–147 (2020)
  * [6] Fortunato, S.: Community detection in graphs. Physics reports 486(3-5), 75–174 (2010)
  * [7] Fujitsu: Digital annealer. https://www.fujitsu.com/global/services/business-services/digital-annealer/
  * [8] Girvan, M., Newman, M.E.: Community structure in social and biological networks. Proceedings of the national academy of sciences 99(12), 7821–7826 (2002)
  * [9] Guimera, R., Amaral, L.A.N.: Functional cartography of complex metabolic networks. nature 433(7028), 895–900 (2005)
  * [10] Hagberg, A., Swart, P., S Chult, D.: Exploring network structure, dynamics, and function using networkx. Tech. rep., Los Alamos National Lab.(LANL), Los Alamos, NM (United States) (2008)
  * [11] Javed, M.A., Younis, M.S., Latif, S., Qadir, J., Baig, A.: Community detection in networks: A multidisciplinary review. Journal of Network and Computer Applications 108, 87–111 (2018)
  * [12] Johnson, S.C.: Hierarchical clustering schemes. Psychometrika 32(3), 241–254 (1967)
  * [13] Knuth, D.E.: The Stanford GraphBase: a platform for combinatorial computing. AcM Press New York (1993)
  * [14] Krogan, N.J., Cagney, G., Yu, H., Zhong, G., Guo, X., Ignatchenko, A., Li, J., Pu, S., Datta, N., Tikuisis, A.P., et al.: Global landscape of protein complexes in the yeast saccharomyces cerevisiae. Nature 440(7084), 637–643 (2006)
  * [15] Kumar, V., Bass, G., Tomlin, C., Dulny, J.: Quantum annealing for combinatorial clustering. Quantum Information Processing 17(2), 39 (2018)
  * [16] Lancichinetti, A., Fortunato, S., Radicchi, F.: Benchmark graphs for testing community detection algorithms. Physical review E 78(4), 046110 (2008)
  * [17] Leskovec, J., Krevl, A.: SNAP Datasets: Stanford large network dataset collection. http://snap.stanford.edu/data (Jun 2014)
  * [18] Liu, X., Ushijima-Mwesigwa, H., Mandal, A., Upadhyay, S., Safro, I., Roy, A.: On modeling local search with special-purpose combinatorial optimization hardware. arXiv preprint arXiv:1911.09810 (2019)
  * [19] Meredith, G.H.: Regular n-valent n-connected nonhamiltonian non-n-edge-colorable graphs. Journal of Combinatorial Theory, Series B 14(1), 55–60 (1973)
  * [20] Naghsh, Z., Javad-Kalbasi, M., Valaee, S.: Digitally annealed solution for the maximum clique problem with critical application in cellular v2x. In: ICC. pp. 1–7. IEEE (2019)
  * [21] Negre, C.F., Ushijima-Mwesigwa, H., Mniszewski, S.M.: Detecting multiple communities using quantum annealing on the d-wave system. Plos one 15(2), e0227538 (2020)
  * [22] Newman, M.E.: Detecting community structure in networks. The European physical journal B 38(2), 321–330 (2004)
  * [23] Newman, M.E.: Spectral methods for community detection and graph partitioning. Physical Review E 88(4), 042822 (2013)
  * [24] Newman, M.E.: Community detection in networks: Modularity optimization and maximum likelihood are equivalent. arXiv preprint arXiv:1606.02319 (2016)
  * [25] Newman, M.E.: Equivalence between modularity optimization and maximum likelihood methods for community detection. Physical Review E 94(5), 052315 (2016)
  * [26] Newman, M.E., Girvan, M.: Finding and evaluating community structure in networks. Physical review E 69(2), 026113 (2004)
  * [27] Palla, G., Derényi, I., Farkas, I., Vicsek, T.: Uncovering the overlapping community structure of complex networks in nature and society. nature 435(7043), 814–818 (2005)
  * [28] Rahman, M.T., Han, S., Tadayon, N., Valaee, S.: Ising model formulation of outlier rejection, with application in wifi based positioning. In: ICASSP. pp. 4405–4409. IEEE (2019)
  * [29] Reichardt, J., Bornholdt, S.: Detecting fuzzy community structures in complex networks with a potts model. Physical Review Letters 93(21), 218701 (2004)
  * [30] Ruan, Y., Fuhry, D., Parthasarathy, S.: Efficient community detection in large networks using content and links. In: Proceedings of the 22nd international conference on World Wide Web. pp. 1089–1098 (2013)
  * [31] Shaydulin, R., Ushijima-Mwesigwa, H., Safro, I., Mniszewski, S., Alexeev, Y.: Community detection across emerging quantum architectures. arXiv preprint arXiv:1810.07765 (2018)
  * [32] Shaydulin, R., Ushijima-Mwesigwa, H., Safro, I., Mniszewski, S., Alexeev, Y.: Network community detection on small quantum computers. Advanced Quantum Technologies 2(9), 1900029 (2019)
  * [33] Traag, V.A., Waltman, L., van Eck, N.J.: From louvain to leiden: guaranteeing well-connected communities. Scientific reports 9(1), 1–12 (2019)
  * [34] Ushijima-Mwesigwa, H., Negre, C.F., Mniszewski, S.M.: Graph partitioning using quantum annealing on the d-wave system. In: Proceedings of the Second International Workshop on Post Moores Era Supercomputing. pp. 22–29 (2017)
  * [35] Ushijima-Mwesigwa, H., Shaydulin, R., Negre, C.F., Mniszewski, S.M., Alexeev, Y., Safro, I.: Multilevel combinatorial optimization across quantum architectures. arXiv preprint arXiv:1910.09985 (2019)
  * [36] Van Rossum, G., et al.: Python programming language. In: USENIX annual technical conference. vol. 41, p. 36 (2007)
  * [37] Zachary, W.W.: An information flow model for conflict and fission in small groups. Journal of anthropological research 33(4), 452–473 (1977)