Statistics
| Branch: | Revision:

root / fiddle / heuristic-betweenness-centrality / centrality_biconnected.py @ 82c7c698

History | View | Annotate | Download (16.2 KB)

1 8e44cb4f Quynh PX Nguyen
# Implement the section IV in Puzis 2012 paper
2 82c7c698 Quynh PX Nguyen
# straight_lineent the section IV in Puzis 2012 paper
3 8e44cb4f Quynh PX Nguyen
# Heuristic Betweenness Centrality - Partitioning to Bi-connected Components
4
5 181e7c50 Quynh PX Nguyen
import os
6 739fe075 Quynh PX Nguyen
import sys
7 3c6ce57c Quynh PX Nguyen
import pprint
8 181e7c50 Quynh PX Nguyen
import networkx as nx
9
import betweenness_centrality as centrality
10 8e44cb4f Quynh PX Nguyen
import utility
11 181e7c50 Quynh PX Nguyen
from pdb import set_trace as debugger
12
13
MAIN_CODE_DIR = os.environ.get('MAIN_CODE_DIR', '')
14
15
16
class HeuristicBetweennessCentrality():
17
    def __init__(self, subgraphs, bicomponents, cutpoints, num_vertices, link_weight, traffic_matrix):
18
        self.subgraphs = subgraphs
19
        self.bicomponents = bicomponents
20
        self.cutpoints = cutpoints
21
        self.num_vertices = num_vertices
22
        self.link_weight = link_weight
23
        self.traffic_matrix = traffic_matrix
24
25
        self.bc_components = list()
26
        self.calculate_bc_non_cutpoint()
27
        self.calculate_bc_cutpoint()
28
29
        self.bc = dict()
30
        self.finalize()
31
32
    def calculate_bc_non_cutpoint(self):
33 3c6ce57c Quynh PX Nguyen
        """BC for non cutpoint
34 181e7c50 Quynh PX Nguyen
        """
35
        for i, subgraphs in enumerate(self.subgraphs):
36
            traffic_matrix = self.traffic_matrix[i]
37 3c6ce57c Quynh PX Nguyen
            results = centrality.weight_betweenness_centrality(subgraphs, traffic_matrix)
38 181e7c50 Quynh PX Nguyen
            self.bc_components.append(results)
39
40
    def calculate_bc_cutpoint(self):
41
        self.bc_cutpoints = dict()
42
43
        bc_inter = dict()
44
        for v in self.cutpoints:
45
            inter = 0
46
            for i, comp in enumerate(self.bicomponents):
47
                if v in comp:
48
                    inter += self.link_weight.get_link_weight(i, v
49
                        ) * self.link_weight.get_reverse_link_weight(i, v)
50
            bc_inter[v] = inter
51
52
        print 'XXXX bc_components = %s' % self.bc_components
53
        print 'XXXX inter = %s' % bc_inter
54
55
        for v in self.cutpoints:
56
            bc_locally = 0
57
            for i, comp in enumerate(self.bicomponents):
58
                if v in comp:
59
                    bc_locally += self.bc_components[i][v]
60 3c6ce57c Quynh PX Nguyen
            # self.bc_cutpoints[v] = bc_locally - bc_inter[v]
61
            # TODO: do not minus the bc_inter
62
            self.bc_cutpoints[v] = bc_locally
63 181e7c50 Quynh PX Nguyen
64
    def finalize(self):
65
        # Add the bc for non cutpoint vertices
66
        for bc_component in self.bc_components:
67
            for key, value in bc_component.iteritems():
68
                if key not in self.bc:
69
                    self.bc[key] = value
70
71
        # Add the bc for cutpoint vertices
72
        for key, value in self.bc_cutpoints.iteritems():
73
            self.bc[key] = value
74
75 3c6ce57c Quynh PX Nguyen
        print '*** betweenness = %s' % self.bc
76 181e7c50 Quynh PX Nguyen
        # Rescale the bc according to the original graph
77 3c6ce57c Quynh PX Nguyen
        factor = 1.0 / ((self.num_vertices - 1) * (self.num_vertices - 2))
78 82c7c698 Quynh PX Nguyen
        # TODO: check the scaling factor, how much should it be?
79 3c6ce57c Quynh PX Nguyen
        # factor = 2.0 / (self.num_vertices*self.num_vertices - 3 * self.num_vertices + 2)
80 181e7c50 Quynh PX Nguyen
        for key, value in self.bc.iteritems():
81
            self.bc[key] = value * factor
82
83
    def write(self, file_suffix=''):
84
        filepath = '/output/heuristic_%s.csv'  % file_suffix
85
        with open(MAIN_CODE_DIR + filepath, 'w') as output:
86
            for node, centrality in self.bc.iteritems():
87
                output.write('%s, %s\n' % (node, centrality))
88
89
    def __str__(self):
90
        return str(self.bc)
91
92
93
class TrafficMatrix():
94
    def __init__(self, bicomponents, cutpoints, link_weight):
95
        self.bicomponents = bicomponents
96
        self.cutpoints = cutpoints
97
        self.num_components = len(bicomponents)
98
        self.link_weight = link_weight
99
100
        self.h = list()
101
        self.generate_empty_traffic_matrix()
102
        self.generate_traffic_matrix()
103
104
    def generate_empty_traffic_matrix(self):
105
        for i in range(self.num_components):
106
            l = len(self.bicomponents[i])
107
            matrix = [[1 for x in range(l)] for y in range(l)]
108
            # update the main diagonal
109
            for x in range(l):
110
                matrix[x][x] = 0
111
112
            self.h.append(matrix)
113
114
    def generate_traffic_matrix(self):
115
        # Update the value when one vertex is a cut-point, another vertex is not a cut-point
116
        for i, components in enumerate(self.bicomponents):
117
            normal_points = components.difference(self.cutpoints)
118
            cutpoints = self.cutpoints.intersection(components)
119
120
            for cutpoint in cutpoints:
121
                for normal_point in normal_points:
122
                    communication_intensity = self.link_weight.get_reverse_link_weight(i, cutpoint) + 1
123
                    self.update(i, cutpoint, normal_point, communication_intensity)
124
125
        # Update the value when both vertices are cut-points
126
        for i, components in enumerate(self.bicomponents):
127
            cutpoints = list(self.cutpoints.intersection(components))
128
            len_cutpoints = len(cutpoints)
129
            if len_cutpoints > 1:
130
                for k in range(len_cutpoints - 1):
131
                    for l in range(1, len_cutpoints):
132 8e44cb4f Quynh PX Nguyen
                        if k == l:
133
                            continue
134 181e7c50 Quynh PX Nguyen
                        communication_intensity = (
135
                            self.link_weight.get_reverse_link_weight(i, cutpoints[k]) + 1) * (
136
                            self.link_weight.get_reverse_link_weight(i, cutpoints[l]) + 1
137
                        )
138
                        self.update(i, cutpoints[k], cutpoints[l], communication_intensity)
139
140
    def simple_update(self, comp_index, x_pos, y_pos, value):
141
        self.h[comp_index][x_pos][y_pos] = value
142
        # to keep the symmetric property of Traffic Matrix
143
        self.h[comp_index][y_pos][x_pos] = value
144
145
    def update(self, comp_index, x, y, value):
146 82c7c698 Quynh PX Nguyen
        comp = sorted(self.bicomponents[comp_index])
147 181e7c50 Quynh PX Nguyen
        try:
148
            x_pos = list(comp).index(x)
149
            y_pos = list(comp).index(y)
150
        except:
151
            debugger()
152
            a = 2
153
154
        self.simple_update(comp_index, x_pos, y_pos, value)
155
156
    def __str__(self):
157
        return str(self.h)
158
159
    def __getitem__(self, key):
160
        return self.h[key]
161
162
163
class LinkWeight():
164
    def __init__(self, graph, bicomponents, cutpoints):
165
        self.num_vertices = nx.number_of_nodes(graph)
166
        self.bicomponents = bicomponents
167
        self.num_components = len(bicomponents)
168
        self.cutpoints = cutpoints
169
170
        self.Dv_B = [dict() for i in range(self.num_components)] # link weight
171 24a4abf9 Quynh PX Nguyen
        self.compute_component_tree_weight()
172 181e7c50 Quynh PX Nguyen
173
        self.reverse_Dv_B = [dict() for i in range(self.num_components)] # reverse link weight
174
        self.generate_reverse_link_weight()
175
176
    def _components_sharing_cutpoint(self, B_cutpoints, point):
177
        indices = list()
178
        for i, cp in enumerate(B_cutpoints):
179
            if point in cp:
180
                indices.append(i)
181
182
        return indices
183
184
    def get_link_weight(self, comp_index, point):
185
        Dv_B_comp = self.Dv_B[comp_index]
186
187
        if point in Dv_B_comp:
188
            return Dv_B_comp[point]
189
        else:
190
            return 0
191
192 3c6ce57c Quynh PX Nguyen
    def set_link_weight(self, comp_index, point, value):
193
        self.Dv_B[comp_index][point] = value
194
195 181e7c50 Quynh PX Nguyen
    def get_reverse_link_weight(self, comp_index, point):
196
        reverse_Dv_B_comp = self.reverse_Dv_B[comp_index]
197
198
        if point in reverse_Dv_B_comp:
199
            return reverse_Dv_B_comp[point]
200
        else:
201
            return 0
202
203
    def generate_link_weight(self):
204
        # How many cutpoints does this component have
205
        B_plus_v = list()
206
        B_cutpoints = list() # number of cutpoints in component B
207
        for comp in self.bicomponents:
208
            points = comp.intersection(self.cutpoints)
209
            B_plus_v.append(comp.difference(points))
210
            B_cutpoints.append(points)
211
212
        len_B_plus_v = [len(x) for x in B_plus_v]
213
        len_B_cutpoints = [len(x) for x in B_cutpoints]
214
215
        # Calculate the Dv_B
216
217
        # For the leaf in the block-cut tree
218
        level = 1
219
        for (i, cp) in enumerate(B_cutpoints):
220
            if len_B_cutpoints[i] == level:
221
                point = list(cp)[0] # there is only 1 element
222
                self.Dv_B[i][point] = len_B_plus_v[i]
223
224
        # For other nodes in the block-cut tree
225
        level += 1
226
        while level <= max(len_B_cutpoints):
227 3c6ce57c Quynh PX Nguyen
            if level == 3:
228
                debugger()
229 181e7c50 Quynh PX Nguyen
            for (i, cp) in enumerate(B_cutpoints):
230
                if len_B_cutpoints[i] == level:
231
232
                    for point in cp:
233 3c6ce57c Quynh PX Nguyen
                        # 1st way
234 181e7c50 Quynh PX Nguyen
                        shared_comp_indices = self._components_sharing_cutpoint(B_cutpoints, point)
235
                        shared_comp_indices.remove(i)
236
                        weight_shared_comp = list()
237
238
                        for index in shared_comp_indices:
239
                            weight_shared_comp.append(self.get_link_weight(index, point))
240
241 3c6ce57c Quynh PX Nguyen
                        weight = self.num_vertices - 1 - sum(weight_shared_comp)
242 181e7c50 Quynh PX Nguyen
243
                        self.Dv_B[i][point] = weight
244
245
            level += 1
246
247
    def generate_reverse_link_weight(self):
248
        for i, Dv_B_i in enumerate(self.Dv_B):
249
            for key, value in Dv_B_i.iteritems():
250
                self.reverse_Dv_B[i][key] = self.num_vertices - 1 - self.get_link_weight(i, key)
251
252 3c6ce57c Quynh PX Nguyen
    def compute_component_tree_weight(self):
253
        """Follows exactly the Algorithm 1 [Puzis 2012]
254
        """
255
        B_cutpoints = list() # number of cutpoints in component B
256
        for comp in self.bicomponents:
257
            points = comp.intersection(self.cutpoints)
258
            B_cutpoints.append(points)
259
260
        Q = self._inititalize_component_tree_weight()
261
        while Q:
262
            pair = Q.pop(0)
263
            if pair['type'] == 'component_vertex_pair':
264
                B_index = pair['value'][0]
265
                v = pair['value'][1]
266
                size = len(self.bicomponents[B_index]) - 1;
267
                all_cutpoints = self.bicomponents[B_index].intersection(self.cutpoints)
268
                all_cutpoints.remove(v)
269
270
                for cp in all_cutpoints:
271
                    if self.get_link_weight(B_index, cp) != -1:
272 24a4abf9 Quynh PX Nguyen
                        size += self.num_vertices - self.get_link_weight(B_index, cp) - 1
273 3c6ce57c Quynh PX Nguyen
274 739fe075 Quynh PX Nguyen
                link_weight = size
275
                self._verify_link_weight(B_index, v, link_weight)
276
                self.set_link_weight(B_index, v, link_weight)
277 3c6ce57c Quynh PX Nguyen
278
                # update Q
279
                Q = self._find_unknown_weight_wrt_cutpoint(v, Q)
280
281
            if pair['type'] == 'vertex_component_pair':
282 24a4abf9 Quynh PX Nguyen
                size = 0
283 3c6ce57c Quynh PX Nguyen
                B_index = pair['value'][0]
284
                v = pair['value'][1]
285
                shared_comp_indices = self._components_sharing_cutpoint(B_cutpoints, v)
286
                shared_comp_indices.remove(B_index)
287
288
289
                for i in shared_comp_indices:
290
                    if self.get_link_weight(i, v) != - 1:
291
                        size += self.get_link_weight(i, v)
292
293 739fe075 Quynh PX Nguyen
                link_weight = self.num_vertices - 1 - size
294
                self._verify_link_weight(B_index, v, link_weight)
295
                self.set_link_weight(B_index, v, link_weight)
296 3c6ce57c Quynh PX Nguyen
297
                # update Q
298
                Q = self._find_unknown_weight_wrt_component(B_index, Q)
299
300 739fe075 Quynh PX Nguyen
    def _verify_link_weight(self, B_index, v, value):
301
        """ If the old_value exist in self.Dv_B, then it must be equal to new value
302

303
        Otherwise, do nothing
304
        """
305
        old_value = self.get_link_weight(B_index, v)
306 82c7c698 Quynh PX Nguyen
307
        if old_value != -1: # -1 is unknown
308 739fe075 Quynh PX Nguyen
            if old_value != value:
309
                print "BUGS FOUND in _verify_link_weight()"
310
                sys.exit()
311
312 24a4abf9 Quynh PX Nguyen
313 3c6ce57c Quynh PX Nguyen
    def _inititalize_component_tree_weight(self):
314
        Q = []
315
        for i, comp in enumerate(self.bicomponents):
316
            current_cutpoints = self.cutpoints.intersection(comp)
317
            if len(current_cutpoints) == 1:
318
                Q.append({
319
                    'type': 'component_vertex_pair',
320
                    'value': (i, list(current_cutpoints)[0]) # (B_i, cutpoint) = (i-th component, the cutpoint name)
321
                    })
322
            for cp in current_cutpoints:
323
                self.set_link_weight(i, cp, -1)
324
        return Q
325
326
    def _find_unknown_weight_wrt_cutpoint(self, cutpoint, Q):
327 24a4abf9 Quynh PX Nguyen
        # Cut-point v such that the weights of all but one of its links in T are already computed.
328
        num_of_uncomputed_weight = 0
329
        uncomputed_component_index = []
330
331 3c6ce57c Quynh PX Nguyen
        for i, Dv_B_comp in enumerate(self.Dv_B):
332 24a4abf9 Quynh PX Nguyen
            if cutpoint in Dv_B_comp:
333
                if Dv_B_comp[cutpoint] == -1:
334
                    num_of_uncomputed_weight += 1
335
                    uncomputed_component_index.append(i)
336
        if num_of_uncomputed_weight == 1:
337
            pair = {
338
                'type': 'vertex_component_pair',
339
                'value': (uncomputed_component_index.pop(), cutpoint)
340
            }
341
            Q.append(pair)
342
        return Q
343
344
    def _find_unknown_weight_wrt_component(self, comp_index, Q):
345
        # Component B such that weights of all but one of its links in T are already computed.
346
        Dv_B_comp = self.Dv_B[comp_index]
347
        values = Dv_B_comp.values()
348
349
        # Check if -1 value appear only 1 time
350
        flag = False
351
        minus_one_value = [x for x in values if x == -1]
352
        if len(minus_one_value) == 1:
353 3c6ce57c Quynh PX Nguyen
            for cp, value in Dv_B_comp.iteritems():
354
                if value == -1:
355
                    pair = {
356 24a4abf9 Quynh PX Nguyen
                        'type': 'component_vertex_pair',
357
                        'value': (comp_index, cp)
358 3c6ce57c Quynh PX Nguyen
                    }
359
                    Q.append(pair)
360
        return Q
361
362 181e7c50 Quynh PX Nguyen
    def __str__(self):
363
        return str(self.Dv_B)
364
365
def analyse_biconnected_component(graph):
366
    bicomponents = list(nx.biconnected_components(graph))
367
    print bicomponents
368
    print [len(x) for x in bicomponents]
369
370 3c6ce57c Quynh PX Nguyen
def find_heuristic_bc(filepath):
371 82c7c698 Quynh PX Nguyen
    filename = os.path.splitext(os.path.basename(filepath))[0]
372 739fe075 Quynh PX Nguyen
373 3c6ce57c Quynh PX Nguyen
    pp = pprint.PrettyPrinter(indent=4)
374 181e7c50 Quynh PX Nguyen
375
    graph = nx.read_weighted_edgelist(filepath)
376
377 739fe075 Quynh PX Nguyen
    if not nx.is_connected(graph):
378
        print "Graph is not connected"
379
        # sys.exit()
380
    else:
381
        print "Graph is connected"
382
383 181e7c50 Quynh PX Nguyen
    # Find biconnected components
384
    analyse_biconnected_component(graph)
385
386
    subgraphs = list(nx.biconnected_component_subgraphs(graph))
387
    bicomponents = list(nx.biconnected_components(graph))
388
    component_edges = list(nx.biconnected_component_edges(graph))
389
    cutpoints = set(nx.articulation_points(graph))
390
391
    num_vertices = nx.number_of_nodes(graph)
392
393
    lw = LinkWeight(graph, bicomponents, cutpoints)
394 3c6ce57c Quynh PX Nguyen
395 24a4abf9 Quynh PX Nguyen
    tm = TrafficMatrix(bicomponents, cutpoints, lw)
396 181e7c50 Quynh PX Nguyen
397 24a4abf9 Quynh PX Nguyen
    bc = HeuristicBetweennessCentrality(subgraphs, bicomponents, cutpoints, num_vertices, lw, tm)
398
    bc.write(file_suffix)
399 739fe075 Quynh PX Nguyen
    print bc
400 181e7c50 Quynh PX Nguyen
401 82c7c698 Quynh PX Nguyen
def test_isomorphic_graphs(filepath1, filepath2):
402
    graph1 = nx.read_weighted_edgelist(filepath1)
403
    graph2 = nx.read_weighted_edgelist(filepath2)
404
    print "Isomorphic = %s" % nx.is_isomorphic(graph1, graph2)
405
406
407 3c6ce57c Quynh PX Nguyen
if __name__ == '__main__':
408 739fe075 Quynh PX Nguyen
    utility.initialize()
409 82c7c698 Quynh PX Nguyen
    f1 = MAIN_CODE_DIR + '/input/9_vertices_green.edges'
410
    f2 = MAIN_CODE_DIR + '/input/9_vertices_blue.edges'
411
    test_isomorphic_graphs(f1, f2)
412 739fe075 Quynh PX Nguyen
413 3c6ce57c Quynh PX Nguyen
    # filepath = MAIN_CODE_DIR + '/input/simple_house.edges'
414
    # filepath = MAIN_CODE_DIR + '/input/simple2.edges'
415
    # filepath = MAIN_CODE_DIR + '/input/simple.edges'
416 739fe075 Quynh PX Nguyen
    # filepath = MAIN_CODE_DIR + '/input/ninux_simplified.edges'
417
    # filepath = MAIN_CODE_DIR + '/input/ninux_simplified_connected.edges'
418 82c7c698 Quynh PX Nguyen
    # filepath = MAIN_CODE_DIR + '/input/ninux_simplified_connected2.edges'
419
    # filepath = MAIN_CODE_DIR + '/input/ninux_simplified_connected3.edges'
420
    # filepath = MAIN_CODE_DIR + '/input/9_vertices_green.edges'
421
    # filepath = MAIN_CODE_DIR + '/input/9_vertices_blue.edges'
422
    # filepath = MAIN_CODE_DIR + '/input/simple_loop.edges'
423
    # filepath = MAIN_CODE_DIR + '/input/straight_line.edges'
424 739fe075 Quynh PX Nguyen
    # filepath = MAIN_CODE_DIR + '/input/ninux_connected.edges'
425
    # filepath = MAIN_CODE_DIR + '/input/ninux_unweighted.edges'
426
    # filepath = MAIN_CODE_DIR + '/input/ninux_unweighted_simplified.edges'
427
    # filepath = MAIN_CODE_DIR + '/input/ninux_unweighted_simplified_connected.edges'
428 3c6ce57c Quynh PX Nguyen
    # filepath = MAIN_CODE_DIR + '/input/simple3.edges'
429 24a4abf9 Quynh PX Nguyen
    # filepath = MAIN_CODE_DIR + '/input/simple4.edges'
430 739fe075 Quynh PX Nguyen
    # filepath = MAIN_CODE_DIR + '/input/simple5.edges'
431 82c7c698 Quynh PX Nguyen
432
    filepath = MAIN_CODE_DIR + '/input/ninux.edges'
433 3c6ce57c Quynh PX Nguyen
    file_suffix = 'edge_list'
434
435 82c7c698 Quynh PX Nguyen
    # Weight betweenness centrality
436
    utility.weight_betweenness_centrality(filepath, file_suffix)
437 739fe075 Quynh PX Nguyen
438 82c7c698 Quynh PX Nguyen
    # Brandess betweenness centrality
439
    utility.brandes_betweenness_centrality(filepath, file_suffix)
440 3c6ce57c Quynh PX Nguyen
441 82c7c698 Quynh PX Nguyen
    # Heuristic BC
442
    find_heuristic_bc(filepath)