Mercurial > repos > shellac > sam_consensus_v3
comparison env/lib/python3.9/site-packages/networkx/algorithms/distance_measures.py @ 0:4f3585e2f14b draft default tip
"planemo upload commit 60cee0fc7c0cda8592644e1aad72851dec82c959"
author | shellac |
---|---|
date | Mon, 22 Mar 2021 18:12:50 +0000 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:4f3585e2f14b |
---|---|
1 """Graph diameter, radius, eccentricity and other properties.""" | |
2 | |
3 import networkx as nx | |
4 from networkx.utils import not_implemented_for | |
5 | |
6 __all__ = [ | |
7 "extrema_bounding", | |
8 "eccentricity", | |
9 "diameter", | |
10 "radius", | |
11 "periphery", | |
12 "center", | |
13 "barycenter", | |
14 "resistance_distance", | |
15 ] | |
16 | |
17 | |
18 def extrema_bounding(G, compute="diameter"): | |
19 """Compute requested extreme distance metric of undirected graph G | |
20 | |
21 Computation is based on smart lower and upper bounds, and in practice | |
22 linear in the number of nodes, rather than quadratic (except for some | |
23 border cases such as complete graphs or circle shaped graphs). | |
24 | |
25 Parameters | |
26 ---------- | |
27 G : NetworkX graph | |
28 An undirected graph | |
29 | |
30 compute : string denoting the requesting metric | |
31 "diameter" for the maximal eccentricity value, | |
32 "radius" for the minimal eccentricity value, | |
33 "periphery" for the set of nodes with eccentricity equal to the diameter | |
34 "center" for the set of nodes with eccentricity equal to the radius | |
35 | |
36 Returns | |
37 ------- | |
38 value : value of the requested metric | |
39 int for "diameter" and "radius" or | |
40 list of nodes for "center" and "periphery" | |
41 | |
42 Raises | |
43 ------ | |
44 NetworkXError | |
45 If the graph consists of multiple components | |
46 | |
47 Notes | |
48 ----- | |
49 This algorithm was proposed in the following papers: | |
50 | |
51 F.W. Takes and W.A. Kosters, Determining the Diameter of Small World | |
52 Networks, in Proceedings of the 20th ACM International Conference on | |
53 Information and Knowledge Management (CIKM 2011), pp. 1191-1196, 2011. | |
54 doi: https://doi.org/10.1145/2063576.2063748 | |
55 | |
56 F.W. Takes and W.A. Kosters, Computing the Eccentricity Distribution of | |
57 Large Graphs, Algorithms 6(1): 100-118, 2013. | |
58 doi: https://doi.org/10.3390/a6010100 | |
59 | |
60 M. Borassi, P. Crescenzi, M. Habib, W.A. Kosters, A. Marino and F.W. Takes, | |
61 Fast Graph Diameter and Radius BFS-Based Computation in (Weakly Connected) | |
62 Real-World Graphs, Theoretical Computer Science 586: 59-80, 2015. | |
63 doi: https://doi.org/10.1016/j.tcs.2015.02.033 | |
64 """ | |
65 | |
66 # init variables | |
67 degrees = dict(G.degree()) # start with the highest degree node | |
68 minlowernode = max(degrees, key=degrees.get) | |
69 N = len(degrees) # number of nodes | |
70 # alternate between smallest lower and largest upper bound | |
71 high = False | |
72 # status variables | |
73 ecc_lower = dict.fromkeys(G, 0) | |
74 ecc_upper = dict.fromkeys(G, N) | |
75 candidates = set(G) | |
76 | |
77 # (re)set bound extremes | |
78 minlower = N | |
79 maxlower = 0 | |
80 minupper = N | |
81 maxupper = 0 | |
82 | |
83 # repeat the following until there are no more candidates | |
84 while candidates: | |
85 if high: | |
86 current = maxuppernode # select node with largest upper bound | |
87 else: | |
88 current = minlowernode # select node with smallest lower bound | |
89 high = not high | |
90 | |
91 # get distances from/to current node and derive eccentricity | |
92 dist = dict(nx.single_source_shortest_path_length(G, current)) | |
93 if len(dist) != N: | |
94 msg = "Cannot compute metric because graph is not connected." | |
95 raise nx.NetworkXError(msg) | |
96 current_ecc = max(dist.values()) | |
97 | |
98 # print status update | |
99 # print ("ecc of " + str(current) + " (" + str(ecc_lower[current]) + "/" | |
100 # + str(ecc_upper[current]) + ", deg: " + str(dist[current]) + ") is " | |
101 # + str(current_ecc)) | |
102 # print(ecc_upper) | |
103 | |
104 # (re)set bound extremes | |
105 maxuppernode = None | |
106 minlowernode = None | |
107 | |
108 # update node bounds | |
109 for i in candidates: | |
110 # update eccentricity bounds | |
111 d = dist[i] | |
112 ecc_lower[i] = low = max(ecc_lower[i], max(d, (current_ecc - d))) | |
113 ecc_upper[i] = upp = min(ecc_upper[i], current_ecc + d) | |
114 | |
115 # update min/max values of lower and upper bounds | |
116 minlower = min(ecc_lower[i], minlower) | |
117 maxlower = max(ecc_lower[i], maxlower) | |
118 minupper = min(ecc_upper[i], minupper) | |
119 maxupper = max(ecc_upper[i], maxupper) | |
120 | |
121 # update candidate set | |
122 if compute == "diameter": | |
123 ruled_out = { | |
124 i | |
125 for i in candidates | |
126 if ecc_upper[i] <= maxlower and 2 * ecc_lower[i] >= maxupper | |
127 } | |
128 | |
129 elif compute == "radius": | |
130 ruled_out = { | |
131 i | |
132 for i in candidates | |
133 if ecc_lower[i] >= minupper and ecc_upper[i] + 1 <= 2 * minlower | |
134 } | |
135 | |
136 elif compute == "periphery": | |
137 ruled_out = { | |
138 i | |
139 for i in candidates | |
140 if ecc_upper[i] < maxlower | |
141 and (maxlower == maxupper or ecc_lower[i] > maxupper) | |
142 } | |
143 | |
144 elif compute == "center": | |
145 ruled_out = { | |
146 i | |
147 for i in candidates | |
148 if ecc_lower[i] > minupper | |
149 and (minlower == minupper or ecc_upper[i] + 1 < 2 * minlower) | |
150 } | |
151 | |
152 elif compute == "eccentricities": | |
153 ruled_out = {} | |
154 | |
155 ruled_out.update(i for i in candidates if ecc_lower[i] == ecc_upper[i]) | |
156 candidates -= ruled_out | |
157 | |
158 # for i in ruled_out: | |
159 # print("removing %g: ecc_u: %g maxl: %g ecc_l: %g maxu: %g"% | |
160 # (i,ecc_upper[i],maxlower,ecc_lower[i],maxupper)) | |
161 # print("node %g: ecc_u: %g maxl: %g ecc_l: %g maxu: %g"% | |
162 # (4,ecc_upper[4],maxlower,ecc_lower[4],maxupper)) | |
163 # print("NODE 4: %g"%(ecc_upper[4] <= maxlower)) | |
164 # print("NODE 4: %g"%(2 * ecc_lower[4] >= maxupper)) | |
165 # print("NODE 4: %g"%(ecc_upper[4] <= maxlower | |
166 # and 2 * ecc_lower[4] >= maxupper)) | |
167 | |
168 # updating maxuppernode and minlowernode for selection in next round | |
169 for i in candidates: | |
170 if ( | |
171 minlowernode is None | |
172 or ( | |
173 ecc_lower[i] == ecc_lower[minlowernode] | |
174 and degrees[i] > degrees[minlowernode] | |
175 ) | |
176 or (ecc_lower[i] < ecc_lower[minlowernode]) | |
177 ): | |
178 minlowernode = i | |
179 | |
180 if ( | |
181 maxuppernode is None | |
182 or ( | |
183 ecc_upper[i] == ecc_upper[maxuppernode] | |
184 and degrees[i] > degrees[maxuppernode] | |
185 ) | |
186 or (ecc_upper[i] > ecc_upper[maxuppernode]) | |
187 ): | |
188 maxuppernode = i | |
189 | |
190 # print status update | |
191 # print (" min=" + str(minlower) + "/" + str(minupper) + | |
192 # " max=" + str(maxlower) + "/" + str(maxupper) + | |
193 # " candidates: " + str(len(candidates))) | |
194 # print("cand:",candidates) | |
195 # print("ecc_l",ecc_lower) | |
196 # print("ecc_u",ecc_upper) | |
197 # wait = input("press Enter to continue") | |
198 | |
199 # return the correct value of the requested metric | |
200 if compute == "diameter": | |
201 return maxlower | |
202 elif compute == "radius": | |
203 return minupper | |
204 elif compute == "periphery": | |
205 p = [v for v in G if ecc_lower[v] == maxlower] | |
206 return p | |
207 elif compute == "center": | |
208 c = [v for v in G if ecc_upper[v] == minupper] | |
209 return c | |
210 elif compute == "eccentricities": | |
211 return ecc_lower | |
212 return None | |
213 | |
214 | |
215 def eccentricity(G, v=None, sp=None): | |
216 """Returns the eccentricity of nodes in G. | |
217 | |
218 The eccentricity of a node v is the maximum distance from v to | |
219 all other nodes in G. | |
220 | |
221 Parameters | |
222 ---------- | |
223 G : NetworkX graph | |
224 A graph | |
225 | |
226 v : node, optional | |
227 Return value of specified node | |
228 | |
229 sp : dict of dicts, optional | |
230 All pairs shortest path lengths as a dictionary of dictionaries | |
231 | |
232 Returns | |
233 ------- | |
234 ecc : dictionary | |
235 A dictionary of eccentricity values keyed by node. | |
236 """ | |
237 # if v is None: # none, use entire graph | |
238 # nodes=G.nodes() | |
239 # elif v in G: # is v a single node | |
240 # nodes=[v] | |
241 # else: # assume v is a container of nodes | |
242 # nodes=v | |
243 order = G.order() | |
244 | |
245 e = {} | |
246 for n in G.nbunch_iter(v): | |
247 if sp is None: | |
248 length = nx.single_source_shortest_path_length(G, n) | |
249 L = len(length) | |
250 else: | |
251 try: | |
252 length = sp[n] | |
253 L = len(length) | |
254 except TypeError as e: | |
255 raise nx.NetworkXError('Format of "sp" is invalid.') from e | |
256 if L != order: | |
257 if G.is_directed(): | |
258 msg = ( | |
259 "Found infinite path length because the digraph is not" | |
260 " strongly connected" | |
261 ) | |
262 else: | |
263 msg = "Found infinite path length because the graph is not" " connected" | |
264 raise nx.NetworkXError(msg) | |
265 | |
266 e[n] = max(length.values()) | |
267 | |
268 if v in G: | |
269 return e[v] # return single value | |
270 else: | |
271 return e | |
272 | |
273 | |
274 def diameter(G, e=None, usebounds=False): | |
275 """Returns the diameter of the graph G. | |
276 | |
277 The diameter is the maximum eccentricity. | |
278 | |
279 Parameters | |
280 ---------- | |
281 G : NetworkX graph | |
282 A graph | |
283 | |
284 e : eccentricity dictionary, optional | |
285 A precomputed dictionary of eccentricities. | |
286 | |
287 Returns | |
288 ------- | |
289 d : integer | |
290 Diameter of graph | |
291 | |
292 See Also | |
293 -------- | |
294 eccentricity | |
295 """ | |
296 if usebounds is True and e is None and not G.is_directed(): | |
297 return extrema_bounding(G, compute="diameter") | |
298 if e is None: | |
299 e = eccentricity(G) | |
300 return max(e.values()) | |
301 | |
302 | |
303 def periphery(G, e=None, usebounds=False): | |
304 """Returns the periphery of the graph G. | |
305 | |
306 The periphery is the set of nodes with eccentricity equal to the diameter. | |
307 | |
308 Parameters | |
309 ---------- | |
310 G : NetworkX graph | |
311 A graph | |
312 | |
313 e : eccentricity dictionary, optional | |
314 A precomputed dictionary of eccentricities. | |
315 | |
316 Returns | |
317 ------- | |
318 p : list | |
319 List of nodes in periphery | |
320 | |
321 See Also | |
322 -------- | |
323 barycenter | |
324 center | |
325 """ | |
326 if usebounds is True and e is None and not G.is_directed(): | |
327 return extrema_bounding(G, compute="periphery") | |
328 if e is None: | |
329 e = eccentricity(G) | |
330 diameter = max(e.values()) | |
331 p = [v for v in e if e[v] == diameter] | |
332 return p | |
333 | |
334 | |
335 def radius(G, e=None, usebounds=False): | |
336 """Returns the radius of the graph G. | |
337 | |
338 The radius is the minimum eccentricity. | |
339 | |
340 Parameters | |
341 ---------- | |
342 G : NetworkX graph | |
343 A graph | |
344 | |
345 e : eccentricity dictionary, optional | |
346 A precomputed dictionary of eccentricities. | |
347 | |
348 Returns | |
349 ------- | |
350 r : integer | |
351 Radius of graph | |
352 """ | |
353 if usebounds is True and e is None and not G.is_directed(): | |
354 return extrema_bounding(G, compute="radius") | |
355 if e is None: | |
356 e = eccentricity(G) | |
357 return min(e.values()) | |
358 | |
359 | |
360 def center(G, e=None, usebounds=False): | |
361 """Returns the center of the graph G. | |
362 | |
363 The center is the set of nodes with eccentricity equal to radius. | |
364 | |
365 Parameters | |
366 ---------- | |
367 G : NetworkX graph | |
368 A graph | |
369 | |
370 e : eccentricity dictionary, optional | |
371 A precomputed dictionary of eccentricities. | |
372 | |
373 Returns | |
374 ------- | |
375 c : list | |
376 List of nodes in center | |
377 | |
378 See Also | |
379 -------- | |
380 barycenter | |
381 periphery | |
382 """ | |
383 if usebounds is True and e is None and not G.is_directed(): | |
384 return extrema_bounding(G, compute="center") | |
385 if e is None: | |
386 e = eccentricity(G) | |
387 radius = min(e.values()) | |
388 p = [v for v in e if e[v] == radius] | |
389 return p | |
390 | |
391 | |
392 def barycenter(G, weight=None, attr=None, sp=None): | |
393 r"""Calculate barycenter of a connected graph, optionally with edge weights. | |
394 | |
395 The :dfn:`barycenter` a | |
396 :func:`connected <networkx.algorithms.components.is_connected>` graph | |
397 :math:`G` is the subgraph induced by the set of its nodes :math:`v` | |
398 minimizing the objective function | |
399 | |
400 .. math:: | |
401 | |
402 \sum_{u \in V(G)} d_G(u, v), | |
403 | |
404 where :math:`d_G` is the (possibly weighted) :func:`path length | |
405 <networkx.algorithms.shortest_paths.generic.shortest_path_length>`. | |
406 The barycenter is also called the :dfn:`median`. See [West01]_, p. 78. | |
407 | |
408 Parameters | |
409 ---------- | |
410 G : :class:`networkx.Graph` | |
411 The connected graph :math:`G`. | |
412 weight : :class:`str`, optional | |
413 Passed through to | |
414 :func:`~networkx.algorithms.shortest_paths.generic.shortest_path_length`. | |
415 attr : :class:`str`, optional | |
416 If given, write the value of the objective function to each node's | |
417 `attr` attribute. Otherwise do not store the value. | |
418 sp : dict of dicts, optional | |
419 All pairs shortest path lengths as a dictionary of dictionaries | |
420 | |
421 Returns | |
422 ------- | |
423 list | |
424 Nodes of `G` that induce the barycenter of `G`. | |
425 | |
426 Raises | |
427 ------ | |
428 NetworkXNoPath | |
429 If `G` is disconnected. `G` may appear disconnected to | |
430 :func:`barycenter` if `sp` is given but is missing shortest path | |
431 lengths for any pairs. | |
432 ValueError | |
433 If `sp` and `weight` are both given. | |
434 | |
435 See Also | |
436 -------- | |
437 center | |
438 periphery | |
439 """ | |
440 if sp is None: | |
441 sp = nx.shortest_path_length(G, weight=weight) | |
442 else: | |
443 sp = sp.items() | |
444 if weight is not None: | |
445 raise ValueError("Cannot use both sp, weight arguments together") | |
446 smallest, barycenter_vertices, n = float("inf"), [], len(G) | |
447 for v, dists in sp: | |
448 if len(dists) < n: | |
449 raise nx.NetworkXNoPath( | |
450 f"Input graph {G} is disconnected, so every induced subgraph " | |
451 "has infinite barycentricity." | |
452 ) | |
453 barycentricity = sum(dists.values()) | |
454 if attr is not None: | |
455 G.nodes[v][attr] = barycentricity | |
456 if barycentricity < smallest: | |
457 smallest = barycentricity | |
458 barycenter_vertices = [v] | |
459 elif barycentricity == smallest: | |
460 barycenter_vertices.append(v) | |
461 return barycenter_vertices | |
462 | |
463 | |
464 def _laplacian_submatrix(node, mat, node_list): | |
465 """Removes row/col from a sparse matrix and returns the submatrix | |
466 """ | |
467 j = node_list.index(node) | |
468 n = list(range(len(node_list))) | |
469 n.pop(j) | |
470 | |
471 if mat.shape[0] != mat.shape[1]: | |
472 raise nx.NetworkXError("Matrix must be square") | |
473 elif len(node_list) != mat.shape[0]: | |
474 msg = "Node list length does not match matrix dimentions" | |
475 raise nx.NetworkXError(msg) | |
476 | |
477 mat = mat.tocsr() | |
478 mat = mat[n, :] | |
479 | |
480 mat = mat.tocsc() | |
481 mat = mat[:, n] | |
482 | |
483 node_list.pop(j) | |
484 | |
485 return mat, node_list | |
486 | |
487 | |
488 def _count_lu_permutations(perm_array): | |
489 """Counts the number of permutations in SuperLU perm_c or perm_r | |
490 """ | |
491 perm_cnt = 0 | |
492 arr = perm_array.tolist() | |
493 for i in range(len(arr)): | |
494 if i != arr[i]: | |
495 perm_cnt += 1 | |
496 n = arr.index(i) | |
497 arr[n] = arr[i] | |
498 arr[i] = i | |
499 | |
500 return perm_cnt | |
501 | |
502 | |
503 @not_implemented_for("directed") | |
504 def resistance_distance(G, nodeA, nodeB, weight=None, invert_weight=True): | |
505 """Returns the resistance distance between node A and node B on graph G. | |
506 | |
507 The resistance distance between two nodes of a graph is akin to treating | |
508 the graph as a grid of resistorses with a resistance equal to the provided | |
509 weight. | |
510 | |
511 If weight is not provided, then a weight of 1 is used for all edges. | |
512 | |
513 Parameters | |
514 ---------- | |
515 G : NetworkX graph | |
516 A graph | |
517 | |
518 nodeA : node | |
519 A node within graph G. | |
520 | |
521 nodeB : node | |
522 A node within graph G, exclusive of Node A. | |
523 | |
524 weight : string or None, optional (default=None) | |
525 The edge data key used to compute the resistance distance. | |
526 If None, then each edge has weight 1. | |
527 | |
528 invert_weight : boolean (default=True) | |
529 Proper calculation of resistance distance requires building the | |
530 Laplacian matrix with the reciprocal of the weight. Not required | |
531 if the weight is already inverted. Weight cannot be zero. | |
532 | |
533 Returns | |
534 ------- | |
535 rd : float | |
536 Value of effective resistance distance | |
537 | |
538 Notes | |
539 ----- | |
540 Overview discussion: | |
541 * https://en.wikipedia.org/wiki/Resistance_distance | |
542 * http://mathworld.wolfram.com/ResistanceDistance.html | |
543 | |
544 Additional details: | |
545 Vaya Sapobi Samui Vos, “Methods for determining the effective resistance,” M.S., | |
546 Mathematisch Instituut, Universiteit Leiden, Leiden, Netherlands, 2016 | |
547 Available: `Link to thesis <https://www.universiteitleiden.nl/binaries/content/assets/science/mi/scripties/master/vos_vaya_master.pdf>`_ | |
548 """ | |
549 import numpy as np | |
550 import scipy.sparse | |
551 | |
552 if not nx.is_connected(G): | |
553 msg = "Graph G must be strongly connected." | |
554 raise nx.NetworkXError(msg) | |
555 elif nodeA not in G: | |
556 msg = "Node A is not in graph G." | |
557 raise nx.NetworkXError(msg) | |
558 elif nodeB not in G: | |
559 msg = "Node B is not in graph G." | |
560 raise nx.NetworkXError(msg) | |
561 elif nodeA == nodeB: | |
562 msg = "Node A and Node B cannot be the same." | |
563 raise nx.NetworkXError(msg) | |
564 | |
565 G = G.copy() | |
566 node_list = list(G) | |
567 | |
568 if invert_weight and weight is not None: | |
569 if G.is_multigraph(): | |
570 for (u, v, k, d) in G.edges(keys=True, data=True): | |
571 d[weight] = 1 / d[weight] | |
572 else: | |
573 for (u, v, d) in G.edges(data=True): | |
574 d[weight] = 1 / d[weight] | |
575 # Replace with collapsing topology or approximated zero? | |
576 | |
577 # Using determinants to compute the effective resistance is more memory | |
578 # efficent than directly calculating the psuedo-inverse | |
579 L = nx.laplacian_matrix(G, node_list, weight=weight) | |
580 | |
581 Lsub_a, node_list_a = _laplacian_submatrix(nodeA, L.copy(), node_list[:]) | |
582 Lsub_ab, node_list_ab = _laplacian_submatrix(nodeB, Lsub_a.copy(), node_list_a[:]) | |
583 | |
584 # Factorize Laplacian submatrixes and extract diagonals | |
585 # Order the diagonals to minimize the likelihood over overflows | |
586 # during computing the determinant | |
587 lu_a = scipy.sparse.linalg.splu(Lsub_a, options=dict(SymmetricMode=True)) | |
588 LdiagA = lu_a.U.diagonal() | |
589 LdiagA_s = np.product(np.sign(LdiagA)) * np.product(lu_a.L.diagonal()) | |
590 LdiagA_s *= (-1) ** _count_lu_permutations(lu_a.perm_r) | |
591 LdiagA_s *= (-1) ** _count_lu_permutations(lu_a.perm_c) | |
592 LdiagA = np.absolute(LdiagA) | |
593 LdiagA = np.sort(LdiagA) | |
594 | |
595 lu_ab = scipy.sparse.linalg.splu(Lsub_ab, options=dict(SymmetricMode=True)) | |
596 LdiagAB = lu_ab.U.diagonal() | |
597 LdiagAB_s = np.product(np.sign(LdiagAB)) * np.product(lu_ab.L.diagonal()) | |
598 LdiagAB_s *= (-1) ** _count_lu_permutations(lu_ab.perm_r) | |
599 LdiagAB_s *= (-1) ** _count_lu_permutations(lu_ab.perm_c) | |
600 LdiagAB = np.absolute(LdiagAB) | |
601 LdiagAB = np.sort(LdiagAB) | |
602 | |
603 # Calculate the ratio of determinant, rd = det(Lsub_ab)/det(Lsub_a) | |
604 Ldet = np.product(np.divide(np.append(LdiagAB, [1]), LdiagA)) | |
605 rd = Ldet * LdiagAB_s / LdiagA_s | |
606 | |
607 return rd |