复杂网络作业六:Louvain社区发现算法原理,细节以及实现


前言

这个第五题本身并不难,只是我个人对这个Louvain的算法比较感兴趣。所以,就花的一周时间。可能是因为这是一篇算法型的论文吧。所以,复现难度不算太大。但是,如果不参考网上的一些已经写完的代码其实也会漏掉很多细节。包括现在也是并不确定我写的是否是一定正确。每次更新之后也会在博客里同步更新。如果是你的目标是想要实现那就从一开始慢慢看,在这里有很多实现细节的说明和我个人的理解。如果要代码就直接跳到最后 ,其实我也不能保证一定正确。毕竟这不是ACM的算法题能有一个准确的对错。


一、Louvain是什么?

Louvain是一个用于社区发现的传统算法。这个算法出现于2008年,那么什么是社区发现呢?举个例子:假设我们有一个图,在这个图中的节点是某个镇的所有的人,图中的每一条边代表的是两个人之间的说话的数量。(现实生活中这个可能难以统计)而社区发现就是在这个图中我们需要确定哪些人之间是一个团体,这个所谓的团体可能是同一个小学,也可能是同一个家庭或者是同一个小区等等。论文的原文全称是:Fast unfolding of communities in large networks 作者:Vincent D. Blondel, Jean-Loup Guillaume, Renaud Lambiotte and Etienne Lefebvre

二、算法思路

1.社区划分的合理性

作为一个认真的人,在放出公式之前先说明。论文中给出的用于衡量社区划分合理性的公式并非是最好的。毕竟,十年都过去了,应该有其他的一些公式吧(我猜的)。但是,我能确定的是最好的一种划分方式一定是真实结果(这句话看起来像句费话)。公式如下(其中的几个乘号是我自己加上去的,我想应该不会错。):
Q = 1 2 ∗ m ∑ i , j [ A i j − k i ∗ k j 2 ∗ m ] δ ( c i , c j ) Q=\frac{1}{2*m}\sum_{i,j} [A_{ij} - \frac{k_i*k_j}{2*m}]\delta(c_i,c_j) Q=2m1i,j[Aij2mkikj]δ(ci,cj)
m m m:如果是无向图,图中所有边的权值和。如果是有向图就是图中所有边的权值和的一半。论文中原本的写法是 m = 1 2 ∑ i j A i j m = \frac{1}{2}\sum_{ij}A_{ij} m=21ijAij表面上看是图中所有边的一半。但是对于无向图而言需要把每一条边当作两条有向边来看。
A i j A_{ij} Aij:代表的是节点i和节点j之间的边的权值。
k i k_i ki:是所有指向节点i的边的权值和。注意:就是这一个定义当我们在处理自环的时候就需要注意一个细节。假设有一个节点的图案如下图所示:
在这里插入图片描述
其中权值为4的边是一个自环需要当作两条边来看,就上图的i点而言其中的 k i k_i ki的值是28。
δ ( c i , c j ) \delta(c_i,c_j) δ(ci,cj):就是判断i节点和j节点是否是同一个社区。如果是那么值就是1否则就是0
我这么写可能是有此啰嗦。但是如果你能有耐心看完那么在实现的时候就能避开很多的坑。因为我这个可不是对论文的简单翻译。

2.算法流程

其算法流程一共有两个大步骤,不断的迭代往复。
1.我们先要把每一个单一的节点都看作是一个单独的社区。对于每一个单一的节点,先要把它从它所在的社区中取出来。然后,再尝试依次加入到它相邻的社区中计算其收益 Δ Q \Delta{Q} ΔQ。( Δ Q \Delta{Q} ΔQ的计算方法之后在详细说,本文给出了三种不同的计算方法。)如果,最大的 Δ Q \Delta{Q} ΔQ小于0则将它放回到原先的社区中,否则就加入到收益最大的那一个社区中。对所有的节点都不断的循环反复的进行这一操作。直到所有节点的所处社区都不再发生变化。(不知道有没有人和我一样看原文的时候这一句没看到。)
从这一步的做法中我主观的感觉这个算法的目的是想要最大化Q的。但是,在实际情况中却又不是Q越大就越接近真实划分。详情就要参考论文中的第12篇参考文献了。(我没有看过,详细了解这个公式Q的同学可以说一下)
2. 通过第一步我们已经对图有了一个新的社区划分。于是在这一步中我们会对每一个社区都要把其中的对应的那些节点都缩成一个超节点,而这个超节点就代表了这一个社区中的若干个结点。而边的处理方式与咱们的逻辑十分相符。如果一条边连接的是同一个社区就把它变成对应超节点的一条自环边,如果一条边连接的是两个不同的社区就让这条边连接两个对应的超节点。然后计算一个新的全局Q值。如果这个Q值与之前的Q值相同那么就结束整个算法。(如果是第一次迭代那就不用考虑结束的事了)
为了能够更加形象的理解,我试试看画一个图。
(假装有图,好吧我画不出来。还是会的工具太少。)

3. Δ Q \Delta{Q} ΔQ的计算方式

在这一部分中我会一共会放三种不同的计算方式这三种方式都是基于一个大前提——我们需要尝试将一个节点加入到社区C中。其中第一种是论文中给出的,第二种是其它的一些代码中使用的版,(据说论文的初稿用的也是这个版本)第三种是我根据自己的理解推导得到的。(其正确性有待商榷,推荐使用第二种。)
1. Δ Q = [ ∑ i n + k i , i n 2 ∗ m − ( ∑ t o t + k i 2 ∗ m ) 2 ] − [ ∑ i n 2 ∗ m − ( ∑ t o t 2 ∗ m ) 2 − ( k i 2 ∗ m ) 2 ] \Delta{Q} = [\frac{\sum_{in}+k_{i,in}}{2*m} - (\frac{\sum_{tot} + k_i}{2*m})^2] - [\frac{\sum_{in}}{2*m}-(\frac{\sum_{tot}}{2*m})^2 - (\frac{k_i}{2*m})^2] ΔQ=[2min+ki,in(2mtot+ki)2][2min(2mtot)2(2mki)2]
由于我个人对这个公式的理解并不完全,因此不对这个公式作非常详细的解读。仅仅只对其中的变量做一些简单的解释。
∑ i n \sum_{in} in:社区C的内部的所有边的权值和,注意:如果是无向图则需要每一条边看作是两条有向边,所以,是所有边的权值和的两倍。
k i , i n k_{i,in} ki,in:所有从节点i指向区域C的边的权值和。
m m m:如果是无向图,图中所有边的权值和。如果是有向图就是图中所有边的权值和的一半。
∑ t o t \sum_{tot} tot:所有指向区域C中的节点的边的权值和。刚开始的时候,我把这个理解成了从区域C指向外面的所有边的权值和。因此,踩了一个大坑。(这机种理解的区别就在于对自环的处理不同。)
k i k_i ki:指向节点i的所有边的权值和。(注意对自环的处理。)
在接下了我们就在这个公式的基础上进行划简:
Δ Q = [ ∑ i n 2 ∗ m + k i , i n 2 ∗ m − ( ∑ t o t 2 ∗ m ) 2 − ( 2 ∗ ∑ t o t ∗ k i 4 ∗ m 2 ) − ( k i 2 ∗ m ) 2 ] − [ ∑ i n 2 ∗ m − ( ∑ t o t 2 ∗ m ) 2 − ( k i 2 ∗ m ) 2 ] = k i , i n 2 ∗ m − 2 ∗ ∑ t o t ∗ k i 4 ∗ m 2 = 1 2 ∗ m ∗ ( k i , i n − ∑ t o t ∗ k i m ) \Delta{Q} = [\frac{\sum_{in}} {2*m} + \frac{k_{i,in}}{2*m} - (\frac{\sum_{tot}}{2*m})^2 - (\frac{2*\sum_{tot}*k_i}{4*m^2}) - (\frac{k_i}{2*m})^2] - [\frac{\sum_{in}}{2*m}-(\frac{\sum_{tot}}{2*m})^2 - (\frac{k_i}{2*m})^2] \\ = \frac{k_{i,in}}{2*m} - \frac{2*\sum_{tot}*k_i}{4*m^2} \\ =\frac{1}{2*m}*(k_{i,in} - \frac{\sum_{tot}*k_i}{m}) ΔQ=[2min+2mki,in(2mtot)2(4m22totki)(2mki)2][2min(2mtot)2(2mki)2]=2mki,in4m22totki=2m1(ki,inmtotki)
在上述公式中 1 2 m \frac{1}{2m} 2m1在算法进行的过程中全程是一个常数,不会影响其大小的比较。因此在比较大小的过程中只需要比较 ( k i , i n − ∑ t o t ∗ k i m ) (k_{i,in} - \frac{\sum_{tot}*k_i}{m}) (ki,inmtotki)即可。
2. Δ Q = 2 ∗ k i , i n − ∑ t o t ∗ k i m \Delta{Q}=2*k_{i,in} - \frac{\sum_{tot}*k_i}{m} ΔQ=2ki,inmtotki
这个公式对应的上一个公式就是在 k i , i n k_{i,in} ki,in前多了一个系数2,据说,作者刚发布的第一个版本的论文中用的就是这个后来修改了系数。同时用用使用这个公式作为 Δ Q \Delta{Q} ΔQ的计算分割效果要更好一点。同时在karate的数据集中使用这个方法的效果确实比上一个公式的效果要更好一点。
3.这个是根据我自己的理解推出来的,但是它的正确性我无法验证。我目前只用它在karate的数据集中运行过,不管是最终的Q值比前两个方法的计算得到的Q值要更大一点。而目测效果与第二个相比相差不大。公式如下:
Δ Q = k i , i n 2 ∗ m − ∑ i k i ∗ f ( i , j ) ∗ k j 4 ∗ m 2 \Delta_{Q} = \frac{k_{i,in}}{2*m} - \frac{\sum_{i}{k_{i}*f(i,j)*k_{j}}}{4*m^2} ΔQ=2mki,in4m2ikif(i,j)kj
其中f(i,j) 表示的是i,j之间是否有直接连边。如果有则为1否则为0。

三、代码实现

在给出我自己的实现代码之前,先给大家另一个网站的实现方法。同时这个实现方法也是我在调试的过程中主要参考的一个代码。放在前面是因为推荐大家使用这个,我想这个代码的应该会比我自己实现的更准确吧。
1.code1
code转载自:https://blog.csdn.net/weixin_40308540/article/details/101269508
(使用方法在最后面)

'''
    Implements the Louvain method.
    Input: a weighted undirected graph
    Ouput: a (partition, modularity) pair where modularity is maximum
'''


class PyLouvain:
    '''
        Builds a graph from _path.
        _path: a path to a file containing "node_from node_to" edges (one per line)
    '''

    @classmethod
    def from_file(cls, path):   #从txt中读取一个net
        f = open(path, 'r')
        lines = f.readlines()
        f.close()
        nodes = {}
        edges = []
        for line in lines:
            n = line.split()
            if not n:
                break
            nodes[n[0]] = 1
            nodes[n[1]] = 1
            w = 1
            if len(n) == 3:
                w = int(n[2])
            edges.append(((n[0], n[1]), w))
        # rebuild graph with successive identifiers
        nodes_, edges_ = in_order(nodes, edges)
        print("%d nodes, %d edges" % (len(nodes_), len(edges_)))
        return nodes_, edges_  #此处作了一点修改

    '''
        Builds a graph from _path.
        _path: a path to a file following the Graph Modeling Language specification
    '''

    @classmethod
    def from_gml_file(cls, path):  #从gml中读取一个net
        f = open(path, 'r')
        lines = f.readlines()
        f.close()
        nodes = {}
        edges = []
        current_edge = (-1, -1, 1)
        in_edge = 0
        for line in lines:
            words = line.split()
            if not words:
                break
            if words[0] == 'id':
                nodes[int(words[1])] = 1
            elif words[0] == 'source':
                in_edge = 1
                current_edge = (int(words[1]), current_edge[1], current_edge[2])
            elif words[0] == 'target' and in_edge:
                current_edge = (current_edge[0], int(words[1]), current_edge[2])
            elif words[0] == 'value' and in_edge:
                current_edge = (current_edge[0], current_edge[1], int(words[1]))
            elif words[0] == ']' and in_edge:
                edges.append(((current_edge[0], current_edge[1]), 1))
                current_edge = (-1, -1, 1)
                in_edge = 0
        nodes, edges = in_order(nodes, edges)
        print("%d nodes, %d edges" % (len(nodes), len(edges)))
        return nodes, edges   #此处作了一点修改

    '''
        Initializes the method.
        _nodes: a list of ints
        _edges: a list of ((int, int), weight) pairs
    '''

    def __init__(self, nodes, edges):
        self.nodes = nodes
        self.edges = edges
        # precompute m (sum of the weights of all links in network)
        #            k_i (sum of the weights of the links incident to node i)
        self.m = 0
        self.k_i = [0 for n in nodes]
        self.edges_of_node = {}
        self.w = [0 for n in nodes]
        self.orginer_network = []
        self.orginer_k_i = []
        for e in edges:
            self.m += e[1]
            self.k_i[e[0][0]] += e[1]
            self.k_i[e[0][1]] += e[1]  # there's no self-loop initially
            # save edges by node
            if e[0][0] not in self.edges_of_node:
                self.edges_of_node[e[0][0]] = [e]
            else:
                self.edges_of_node[e[0][0]].append(e)
            if e[0][1] not in self.edges_of_node:
                self.edges_of_node[e[0][1]] = [e]
            elif e[0][0] != e[0][1]:
                self.edges_of_node[e[0][1]].append(e)
        # access community of a node in O(1) time
        self.communities = [n for n in nodes]
        self.actual_partition = []
        self.orginer_k_i = []
        for k in self.k_i:
            self.orginer_k_i.append(k)
    '''
        Applies the Louvain method.
    '''
    def findRoot(self,node):
        for i,community in enumerate(self.actual_partition):
            if(node in community):
                return i
    def my_compute_modularity(self):
        sum_Q = 0
        for edge in self.orginer_network[1]:
            u,v = edge[0]
            w = edge[1]
            if(self.findRoot(u) == self.findRoot(v)):
                sum_Q = (w - (self.orginer_k_i[u]*self.orginer_k_i[v])/(2*self.m))/(2*self.m)
        return sum_Q
    def apply_method(self):
        network = (self.nodes, self.edges)
        self.orginer_network = network
        best_partition = [[node] for node in network[0]]
        best_q = -1
        i = 1
        while 1:
            # print("pass #%d" % i)
            i += 1
            partition = self.first_phase(network)
            #q = self.compute_modularity(partition)
            q = self.my_compute_modularity()
            partition = [c for c in partition if c]
            # print("%s (%.8f)" % (partition, q))
            # clustering initial nodes with partition
            if self.actual_partition:
                actual = []
                for p in partition:
                    part = []
                    for n in p:
                        part.extend(self.actual_partition[n])
                    actual.append(part)
                self.actual_partition = actual
            else:
                self.actual_partition = partition
            print(q)
            if q == best_q:
                break
            network = self.second_phase(network, partition)
            best_partition = partition
            best_q = q
        return (self.actual_partition, best_q)

    '''
        Computes the modularity of the current network.
        _partition: a list of lists of nodes
    '''

    def compute_modularity(self, partition):
        q = 0
        m2 = self.m * 2
        for i in range(len(partition)):
            q += self.s_in[i] / m2 - (self.s_tot[i] / m2) ** 2
        return q

    '''
        Computes the modularity gain of having node in community _c.
        _node: an int
        _c: an int
        _k_i_in: the sum of the weights of the links from _node to nodes in _c
    '''

    def compute_modularity_gain(self, node, c, k_i_in):
        return 2 * k_i_in - self.s_tot[c] * self.k_i[node] / self.m

    '''
        Performs the first phase of the method.
        _network: a (nodes, edges) pair
    '''

    def first_phase(self, network):
        # make initial partition
        best_partition = self.make_initial_partition(network)
        while 1:
            improvement = 0
            for node in network[0]:
                node_community = self.communities[node]
                # default best community is its own
                best_community = node_community
                best_gain = 0
                # remove _node from its community
                best_partition[node_community].remove(node)
                best_shared_links = 0
                for e in self.edges_of_node[node]:
                    if e[0][0] == e[0][1]:
                        continue
                    if e[0][0] == node and self.communities[e[0][1]] == node_community or e[0][1] == node and \
                            self.communities[e[0][0]] == node_community:
                        best_shared_links += e[1]
                self.s_in[node_community] -= 2*(best_shared_links + self.w[node])
                #self.s_tot[node_community] -= self.k_i[node]
                self.s_tot[node_community] -= (self.k_i[node] - 2*best_shared_links)
                self.communities[node] = -1
                communities = {}  # only consider neighbors of different communities
                for neighbor in self.get_neighbors(node):
                    community = self.communities[neighbor]
                    if community in communities:
                        continue
                    communities[community] = 1
                    shared_links = 0
                    for e in self.edges_of_node[node]:
                        if e[0][0] == e[0][1]:
                            continue
                        if e[0][0] == node and self.communities[e[0][1]] == community or e[0][1] == node and \
                                self.communities[e[0][0]] == community:
                            shared_links += e[1]
                    # compute modularity gain obtained by moving _node to the community of _neighbor
                    gain = self.compute_modularity_gain(node, community, shared_links)
                    if gain > best_gain:
                        best_community = community
                        best_gain = gain
                        best_shared_links = shared_links
                # insert _node into the community maximizing the modularity gain
                best_partition[best_community].append(node)
                self.communities[node] = best_community
                self.s_in[best_community] += 2*( best_shared_links + self.w[node])
                #self.s_tot[best_community] += (self.k_i[node])
                self.s_tot[best_community] += (self.k_i[node] - 2*best_shared_links)
                if node_community != best_community:
                    improvement = 1
            if not improvement:
                break
        return best_partition

    '''
        Yields the nodes adjacent to _node.
        _node: an int
    '''

    def get_neighbors(self, node):
        for e in self.edges_of_node[node]:
            if e[0][0] == e[0][1]:  # a node is not neighbor with itself
                continue
            if e[0][0] == node:
                yield e[0][1]
            if e[0][1] == node:
                yield e[0][0]

    '''
        Builds the initial partition from _network.
        _network: a (nodes, edges) pair
    '''

    def make_initial_partition(self, network):
        partition = [[node] for node in network[0]]
        self.s_in = [0 for node in network[0]]
        self.s_tot = [self.k_i[node] for node in network[0]]
        for e in network[1]:
            if e[0][0] == e[0][1]:  # only self-loops
                self.s_in[e[0][0]] += e[1]
                self.s_in[e[0][1]] += e[1]
        return partition

    '''
        Performs the second phase of the method.
        _network: a (nodes, edges) pair
        _partition: a list of lists of nodes
    '''

    def second_phase(self, network, partition):
        nodes_ = [i for i in range(len(partition))]
        # relabelling communities
        communities_ = []
        d = {}
        i = 0
        for community in self.communities:
            if community in d:
                communities_.append(d[community])
            else:
                d[community] = i
                communities_.append(i)
                i += 1
        self.communities = communities_
        # building relabelled edges
        edges_ = {}
        for e in network[1]:
            ci = self.communities[e[0][0]]
            cj = self.communities[e[0][1]]
            try:
                edges_[(ci, cj)] += e[1]
            except KeyError:
                edges_[(ci, cj)] = e[1]
        edges_ = [(k, v) for k, v in edges_.items()]
        # recomputing k_i vector and storing edges by node
        self.k_i = [0 for n in nodes_]
        self.edges_of_node = {}
        self.w = [0 for n in nodes_]
        for e in edges_:
            self.k_i[e[0][0]] += e[1]
            self.k_i[e[0][1]] += e[1]
            if e[0][0] == e[0][1]:
                self.w[e[0][0]] += e[1]
            if e[0][0] not in self.edges_of_node:
                self.edges_of_node[e[0][0]] = [e]
            else:
                self.edges_of_node[e[0][0]].append(e)
            if e[0][1] not in self.edges_of_node:
                self.edges_of_node[e[0][1]] = [e]
            elif e[0][0] != e[0][1]:
                self.edges_of_node[e[0][1]].append(e)
        # resetting communities
        self.communities = [n for n in nodes_]
        return (nodes_, edges_)


'''
    Rebuilds a graph with successive nodes' ids.
    _nodes: a dict of int
    _edges: a list of ((int, int), weight) pairs
'''


def in_order(nodes, edges):
    # rebuild graph with successive identifiers
    nodes = list(nodes.keys())
    nodes.sort()
    i = 0
    nodes_ = []
    d = {}
    for n in nodes:
        nodes_.append(i)
        d[n] = i
        i += 1
    edges_ = []
    for e in edges:
        edges_.append(((d[e[0][0]], d[e[0][1]]), e[1]))
    return (nodes_, edges_)
###################******************---下面是使用方法的一个例子---******************###################
#下面是Louvain的使用方式
import argparse
import networkx as nx
import random
import  matplotlib.pyplot as plt
def getRandomColor():
    '''
    这是随机获取一种颜色。但是,不能保证获取的颜色差距一定很大。
    所以,如果想要直观的看到结果。有时候需要多运行几次。
    '''
    return random.randint(0,255)
def drawNet_Louvain(G,part_list):
    '''
    将社区划分完成的图进行直观的展示。
    '''
    for i in range(len(part_list)):
        for j in range(len(part_list[i])):
            part_list[i][j] += 1
    print(part_list)
    color_list = []
    for part in part_list:
        color = getRandomColor()
        for node in part:
            color_list.append((node,color))
    color_list.sort(key = lambda x: x[0])
    print(color_list)
    color_list = [x[1] for x in color_list]
    plt.figure(figsize=(5, 5))
    print("finish")
    print(len(G.nodes()))
    print(len(color_list))
    pos = nx.spring_layout(G)
    nx.draw(G, with_labels=True, node_color=color_list, pos=pos)
    plt.show()
    plt.savefig(r"filename.png")


def main(option):
    
    data_path = option.data_path#此处需要修改对应的数据文件的路径
    if (data_path[-3:] == "txt"):  #如果文件是在txt中的读取方式
        net_G_nodes,net_G_edges = PyLouvain.from_file(data_path) #使他的方法进行数据读取,返回的是点集和边集
        net_G = nx.read_edgelist(data_path)   #因为他使用的是完全自己实现的代码,无法进行画图展示。所以,需要自己在读入一个networkx的
    elif (data_path[-3:] == "gml"):  #如果文件是在gml中的读取方式
        net_G_nodes,net_G_edges = PyLouvain.from_gml_file(data_path)#使他的方法进行数据读取,返回的是点集和边集
        net_G  = nx.read_gml(data_path,label="id") #因为他使用的是完全自己实现的代码,无法进行画图展示。所以,需要自己在读入一个networkx的
    new_G = PyLouvain(net_G_nodes,net_G_edges)  #使用它的方法构造成一个类,传入的参数依次是点集和边集
    t,d = new_G.apply_method() #应用其中的方法,对社区进行分割。返回值分别是最佳的社区划分,以及对应的Q值。
    drawNet_Louvain(net_G, t) #对分割完成的图进行展示。

if __name__ == "__main__":
    parser = argparse.ArgumentParser()
    parser.add_argument("--data_path", type=str, default="karate.gml",   
                        help="data_path is the path of the net")

    opt = parser.parse_args()
    main(opt)

然后是我根据我自己的理解所实现的算法
2.code

import networkx as nx
import argparse

def CalDeltaQ1(G,node_i,k_i_in,k_i,community_list,community_value,community_tot,
               community_id,neg_community_id,node_edge,edge_list,m,opt):
    '''
    function:用于计算deltaQ
    parameter:
        G:当前的图
        node_i:正在操作的节点的编号
        k_i_in:从i指向区域C的所有边的权值和
        k_i:指向节点i的所有边的权值和(包括自环的边)
        community_list:每一个社区所包含的结点的编号
        community_value:社区内的所有这的权值和
        community_tot:社区指向外部节点的所有边的权值和
        community_id :node_i 所处社区的编号
        neg_community_id:相邻社区的社区编号
        node_edge:node_i这一节点的邻边
        edge_list:图的所有边的列表
        m:全图的所有边的和
        opt:其中有选择deltaQ的计算方式的选项
    return:返回对应的deltaQ
    '''
    sumC = community_value[neg_community_id]
    sum_tot = community_tot[neg_community_id] + community_value[neg_community_id]
    if(opt.calDeltaQfunction == "fun_copy"):  #使用的是参考代码中的公式将其中的2改成1就是优化后的论文中的公式
        return 2*k_i_in - sum_tot*k_i/m
    elif(opt.calDeltaQfunction == "my_fun"):  #使用我自己推导的公式,因为没有验证过它的正确性。所以不推荐使用
        sum_k_i = 0
        delta_Q = (k_i_in)/(2*m)
        k_j = 0
        for edge in G[node_i]:
            if(edge in community_list[neg_community_id]):
                k_j += G.nodes[edge]["node_out"] + G.nodes[edge]["value"]
        delta_Q -= (k_j*k_i)/(4*(m**2))
        return delta_Q
    else:   #原模原样的使用论文中的公式
        m2 = 2*m
        return ((sumC + k_i_in)/m2 - ((sum_tot + k_i)/m2)**2) - (sumC/m2- (sum_tot/m2)**2 - (k_i/m2)**2)

def PutNode2Community(node,community_list,community_value,community_tot,k_i_c,k_i,
                      root_id,m,level):
    '''
    使用并查集的方法进行编号合并
    parameter:
       node:结点 包涵信息
       root_id: 需要把结点加入的社区的id
       m:全图的权值总合
       level:经过的第几次合并(不知道怎么用)
    return:
        无,在过程中会更新node_list中的各个结点的root_list属性值
    '''
    node["root_list"][0] = root_id
    community_list[root_id].append(node["node_id"])
    community_value[root_id] += (node["value"] + 2*k_i_c)
    community_tot[root_id] += (k_i - 2*k_i_c)

def getInitCommunity(G,node_list):  #社区中所包涵的节点,社区value(sumC_in),社区的出度(sumC_tot)
    node2community = []
    node2community_value = []
    node2community_tot = []
    for node in node_list:
        node2community.append([node])
        node2community_value.append(G.nodes[node]["value"])
        node2community_tot.append(G.nodes[node]["node_out"])
    return node2community,node2community_value,node2community_tot
def flag2Merge(G,node_list,m,level,opt):
    '''
    function :对接收的图G进行社区划分,返回划分结果
        node_list:结点列表
        m        :全图的权值和
        level    :迭代的次数
        opt      :里面包含了deltaQ的计算方式
    return:
        社区划分的表示形式。
        "set":每一个社区中所包含的节点的编号。
        "value":社区内部的权值。(社区内的边的两倍)
        "tot":社区内的节点指向,社区外的边的和
    '''
    #对社区的信息进行初始化
    community_list,community_value,community_tot = getInitCommunity(G,node_list)
    while(True):
        promote = False   #是否对于社区划分是否有更新
        for node in node_list:
            #先将node结点移出node所对应的社区  可优化
            community_id = G.nodes[node]["root_list"][0]   #该节点所属的社区的编号
            community_list[community_id].remove(node)      #将节点移出其所属的社区
            k_i_c = 0                                      #从节点指向社区的所有边的权值和
            for edge in G[node]:
                if(edge in community_list[community_id]):
                    k_i_c += G[node][edge]["weight"]

            community_value[community_id] -= (2*k_i_c + G.nodes[node]["value"])  #由于已经将节点移出社区,所以需要更新其区社的value值
            community_tot[community_id] -= (G.nodes[node]["node_out"] - 2*k_i_c) #由于已经将节点移出社区,所以需要更新其区社的tot值
            #在相邻的社区中找最优
            max_delta_Q = 0
            max_k_i_c = k_i_c
            best_community_id = community_id
            #k_i = sum([edge["weight"] for edge in G[node].values()]) #原写法
            k_i = G.nodes[node]["node_out"]  #此处先把k_i 表示为该节点指向其它节点的边的权值和
            #尝试将节点加入到相邻的社区之中
            neg_community_use = set()
            for neg_node in G[node]:
                if(isinstance(neg_node,int)):
                    neg_community = G.nodes[neg_node]["root_list"][0]   #相邻节点所属的社区
                    if(neg_community in neg_community_use):  #如果已经尝试过则不再尝试
                        continue
                    else:
                        neg_community_use.add(neg_community)
                    k_i_c = 0  #节点相相邻节点之间的权值和
                    for edge in G[node]:
                        if(edge in community_list[neg_community]):
                            k_i_c += G[node][edge]["weight"]
                    #计算deltaQ
                    delta_Q = CalDeltaQ1(G,node,k_i_c,k_i + G.nodes[node]["value"],
                                         community_list,community_value,community_tot,
                                         community_id,neg_community,
                                         G[node],
                                         G.edges(),m,opt)
                    #更新最优划分,以及相关信息
                    if(delta_Q > max_delta_Q):
                        max_delta_Q = delta_Q
                        max_k_i_c = k_i_c
                        best_community_id = G.nodes[neg_node]["root_list"][0]
            #将节点加入到最优划分的区间中
            if(max_delta_Q > 0 and best_community_id != community_id):
                promote = True
            #将节点加入到最优的社区中,并且更新其comuniyu的相关属性
            PutNode2Community(G.nodes[node], community_list,community_value,community_tot,max_k_i_c,k_i,
                              best_community_id, m, level)
        if(not promote):
            break

    return {"set":community_list,"value":community_value,"tot":community_tot}
def MergeG2getNewG(G,community):
    '''
    function:根据其社区划分的结果重新生成一个新的图
    G:原先的图结构
    community:社区划分的结果,
        community["set"]——社区中所包含的节点的编号
        community["value"]——社区内部的所有边的权值和的两倍
        community["tot"]——所有社区内部节点指向社区外部节点的边的权值和
    return:返回一个缩点后的图,以及其中对应的边
    '''
    new_G = nx.Graph() #一个新的图
    cnt_community = 0
    #对社区进行遍历
    for i,community_set in enumerate(community["set"]):  #将节点进行合并产生一批新图中的总节点
        #node_id node_set root_list value out
        if(community_set):   #如果在社区中有元素则构造新的节点,加入到图中
            #编号i作为节点编号,community["set"]中的元素就是该社区所包含的对应的元素,将cnt_community作为临时根,community["value"]作为该节点的编号,community["tot"]作为该节点的出度
            #community["value"] + community["tot"]就是指向该社区节点的所有边的权值和
            new_G.add_node(i,node_id = i,node_set = community_set,root_list = [cnt_community],value = community["value"][i],node_out = community["tot"][i])
            cnt_community += 1
    for edge in G.edges(): #将边进行处理加入到适当的新图中
        #weight
        node1,node2 = edge
        root1 = G.nodes[node1]["root_list"][0]
        root2 = G.nodes[node2]["root_list"][0]
        if(root1 == root2): #这条边的一个社区的内部
            #在community["value"]中已经加上,所以不用在这里再加。在PutNode2Community函数中
            #new_G.nodes[root1]["value"] += 2*G[node1][node2]["weight"]
            continue
        else:
            #这条边的两个社区之间,则与其它的边一同进行合并
            if((root1,root2) not in new_G.edges()):
                new_G.add_edge(root1,root2,weight = G[node1][node2]["weight"])
            else:
                new_G[root1][root2]["weight"] += G[node1][node2]["weight"]
    return new_G

def dealSelfLoop(G):
    '''
    对图中的所有的边进行遍历,如果存在自环
    就将其删除,同时一个自环就对那个对应的节点上的值+2
    问题:对于networkx 而言它会自动的删除重边。所以在这一点上对于
    使用networkx对最终的结果会有影响。
    '''

    for node in G:
        cnt_self_node = 0
        flag = False
        for neg_node in G[node]:
            if(node == neg_node): #如果有重边就删除同时更新节点上的值
                flag = True
                cnt_self_node += 2
        if(flag):
            G.remove_edge(node, node)
        G.nodes[node]["value"] = cnt_self_node
    return G

def stateINIT(G):
    for node in G:
        for edge in G[node]:
            G[node][edge]["weight"] = 1   #每一条边都有一个权值,由于是无权图因此权值都需要变为1
    for i,node in enumerate(G.nodes):
        G.nodes[node]["node_id"] = node   #节点的ID
        G.nodes[node]["node_set"] = {node}#社区中所包含的节点(对于第一个图用不上)
        G.nodes[node]["root_list"] = [i]  #每一个节点从属的社区的编号
        G.nodes[node]["node_out"] = sum([G[node][edge]["weight"] for edge in G[node]]) #节点相邻边的权值和

def graphArrange(levelG_list):
    '''
    图的整理,也就是根节点的向下传递
    更新每一层的图的根节点
    '''
    l = -1*len(levelG_list) - 1
    for level in range(-2,l,-1):
        G1 = levelG_list[level + 1]
        G2 = levelG_list[level]
        for node in G1.nodes():
            #获取节点的根,然后将该社区内的所有节点的根也进行更新
            root_id = G1.nodes[node]["root_list"][-1]
        #    print(root_id,":")
            for member_id in G1.nodes[node]["node_set"]:
        #       print(member_id,end=" ")
                G2.nodes[member_id]["root_list"].append(root_id)
        #     print("")
        # print(G1.nodes())
        # print(G2.nodes())

def calGraphQ(G,new_G,m,level,opt):
    '''
    function:计算图在当前的社区划分下的Q值
        G:原图
        new_G:最近的经过缩点后的图
        m:全图的权值和
        level:递归的次数
        opt:其中指定了是使用哪一种Q的计算方法
    return:图的Q值
    '''
    if(opt.calQFunction == "fun_copy"):   #这是参考代码中使用的Q的计算方式
        sum_Q = 0
        for new_node in new_G:
            sum_in = new_G.nodes[new_node]["value"]
            sum_tot = new_G.nodes[new_node]["node_out"]
            sum_Q += (sum_in/(2*m) - (sum_tot/(2*m))**2)
        return sum_Q
    else:  #这是在论文中使用的Q的计算方式
        sum_Q = 1
        for edge in G.edges:
            u,v = edge
            root_u = G.nodes[u]["root_list"][-1]
            root_v = G.nodes[v]["root_list"][-1]
            if(root_u == root_v): #如果u和v同属于一个社区
                #G.nodes[u]["value"]节点内部的边(自环)
                #G.nodes[v]["node_out"] 节点的邻边
                sum_u = G.nodes[u]["value"] + G.nodes[v]["node_out"]
                sum_v = G.nodes[v]["value"] + G.nodes[v]["node_out"]
                sum_Q -= sum_u*sum_v/(4*(m**2))*2   #乘二是因为(u,v)(v,u)要当作两条边来计算
                #sum_Q += (G[u][v]["weight"] - (sum_u)*(sum_v )/(2*m))/(2*m)

        return sum_Q

def Louvain(G,m,opt):
    dealSelfLoop(G)  #对自环进行处理
    stateINIT(G)   #对G添加属性,同时完成初始化
    level = 0   #记录迭代次数
    merageG_list = []  #在这个列表中用于记录在迭代过程中的所有形成的子图
    best_q = -2  #不断更新,获取最终的Q值
    while(True):
        level += 1
        #返回的是对图的一个社区划分,以及每一个社区中的信息
        #其中
        # community["set"]中存储的是每一个社区中所包涵的结点的编号
        # community["value"] 中表示的是每一个社区的权值,其权值主要是自环的2倍
        # sommunity["tot"]中表示的是每一个社区内的结点指向外面的所有边的权值和。(注意:此外不是指向所有社区内的所有边的权值和)
        community = flag2Merge(G,list(G.nodes),m,level,opt)
        #将图的社区划分进行缩点重新组成一个新的图,同时对于新的图添加必要的信息
        new_G = MergeG2getNewG(G,community)
        #将新的图添加到时图列表中
        merageG_list.append(G)
        #将图列表进行整理,只整理根节点的信息
        graphArrange(merageG_list)
        #重新更新使用的图
        G = new_G
        #计算图的Q值
        new_q = calGraphQ(merageG_list[0],G,m,level,opt)
        print(new_q)
        #如果Q值没有变化则,线束算法
        if(new_q == best_q):
            break
        #更新最优的Q值(其实就是最后一种划分的Q值)
        best_q = new_q
    return merageG_list

def CalGraphM1(net_G): #把单个边当作无向图的一条边来处理
    #因为是无权所以边的权值看作是1
    m = len(net_G.edges())
    return m
import random
import matplotlib.pyplot as plt
def getRandomColor():
    return random.randint(0,255)
def drawNet_Louvain(G_Louvain_list,show_level):
    """
    function:将每一次迭代产生的社区划分结果进行展示
    
    """
    show_G = G_Louvain_list[0]
    #print(len(G_Louvain_list[show_level].nodes))
    color_list = []
    color_dict = {}
    for node in show_G.nodes():
        node_color = show_G.nodes[node]["root_list"][show_level]
        if(node_color in color_dict):
            color_list.append(color_dict[node_color])
        else:
            color = getRandomColor()
            color_dict[node_color] = color
            color_list.append(color_dict[node_color])
    plt.figure(figsize=(5, 5))
    #print("finish")
    #print(len(show_G.nodes()))
    #print(len(color_list))
    pos = nx.spring_layout(show_G)
    nx.draw(show_G,with_labels=True,node_color=color_list,pos=pos)
    #plt.show()
    plt.savefig("filename" + str(show_level) + ".png")

def main(option):
    data_path = option.data_path
    if(data_path[-3:] == "txt"):
        net_G = nx.read_edgelist(data_path,
                                 comments='#',
                                 delimiter=None,
                                 create_using=None,
                                 nodetype=int,
                                 data=True,
                                 edgetype=int,
                                 encoding='utf-8',
                                )
    elif(data_path[-3:] == "gml"):
        net_G = nx.read_gml(data_path,
                            label="id")
    m = CalGraphM1(net_G)  #计算全图换边权和
    G_merage_list = Louvain(net_G,m,option)  #进行Louvain算法
    for i in range(len(G_merage_list)):   #画图
        drawNet_Louvain(G_merage_list,i)
if __name__ == "__main__":
    parser = argparse.ArgumentParser()
    parser.add_argument("--data_path",type = str,default="karate.gml",   #数据的路径
                        help="data_path is the path of the net")
    parser.add_argument("--calQFunction",type = str,default="22",   #计算Q的方式
    #Q的选择只能是使用"fun_copy"或者其它,如果是"fun_copy"则使用的是code1中的Q的计算方法。如果是其它则使用的是原论文中的Q的计算方法
                        help="calQFunction is the mode to calQ,should choose from ['fun_copy',other]") 
    parser.add_argument("--calDeltaQfunction",type=str ,default="my_fun",  #计算deltaQ的方式
     #deltaQ的选择只能是使用"fun_copy"或者其它,如果是"fun_copy"则使用的是code1中的deltaQ的计算方法。如果是"my_fun"则使用的是我自己的理解推的一个公式。如果是其它则使用的是原论文中的deltaQ的计算方法
                        help="calDeltaQfunction the mode to calDeltaQ,should choose from ['fun_copy','my_fun',other]")
    opt = parser.parse_args()
    #如果使用的是jupyter需要把上面一行注释掉修改成下面这一行
    #opt = args = parser.parse_args(args=[])
    main(opt)

总结

我已经把我的代码进行过更新了,暂时先这样吧!效果也还算可观。如果,之后有变化那就在更新。另外,好像第五题的第一题就是代码(略)。第二题就是之前的东西(略)。第三题就是跑代码(略)。第四题看情况吧,如果我做了我就一定会发博客。


因为一些原因在此添加数据的下载链接。(上面代码使用的数据) 下载链接:http://www-personal.umich.edu/~mejn/netdata/ 数据名称:Zachary's karate club 数据中的结点较少可以直接看到。

最后,如果有问题可以直接留言咨询。如果我看到了,我一定会回复的

  • 38
    点赞
  • 157
    收藏
    觉得还不错? 一键收藏
  • 51
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 51
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值