Coverage for /pythoncovmergedfiles/medio/medio/usr/local/lib/python3.9/dist-packages/networkx/algorithms/distance_measures.py: 9%
243 statements
« prev ^ index » next coverage.py v7.3.2, created at 2023-10-20 07:00 +0000
1"""Graph diameter, radius, eccentricity and other properties."""
3import networkx as nx
4from networkx.utils import not_implemented_for
6__all__ = [
7 "eccentricity",
8 "diameter",
9 "radius",
10 "periphery",
11 "center",
12 "barycenter",
13 "resistance_distance",
14 "kemeny_constant",
15]
18def _extrema_bounding(G, compute="diameter", weight=None):
19 """Compute requested extreme distance metric of undirected graph G
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).
25 Parameters
26 ----------
27 G : NetworkX graph
28 An undirected graph
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 "eccentricities" for the maximum distance from each node to all other nodes in G
37 weight : string, function, or None
38 If this is a string, then edge weights will be accessed via the
39 edge attribute with this key (that is, the weight of the edge
40 joining `u` to `v` will be ``G.edges[u, v][weight]``). If no
41 such edge attribute exists, the weight of the edge is assumed to
42 be one.
44 If this is a function, the weight of an edge is the value
45 returned by the function. The function must accept exactly three
46 positional arguments: the two endpoints of an edge and the
47 dictionary of edge attributes for that edge. The function must
48 return a number.
50 If this is None, every edge has weight/distance/cost 1.
52 Weights stored as floating point values can lead to small round-off
53 errors in distances. Use integer weights to avoid this.
55 Weights should be positive, since they are distances.
57 Returns
58 -------
59 value : value of the requested metric
60 int for "diameter" and "radius" or
61 list of nodes for "center" and "periphery" or
62 dictionary of eccentricity values keyed by node for "eccentricities"
64 Raises
65 ------
66 NetworkXError
67 If the graph consists of multiple components
68 ValueError
69 If `compute` is not one of "diameter", "radius", "periphery", "center", or "eccentricities".
71 Notes
72 -----
73 This algorithm was proposed in [1]_ and discussed further in [2]_ and [3]_.
75 References
76 ----------
77 .. [1] F. W. Takes, W. A. Kosters,
78 "Determining the diameter of small world networks."
79 Proceedings of the 20th ACM international conference on Information and knowledge management, 2011
80 https://dl.acm.org/doi/abs/10.1145/2063576.2063748
81 .. [2] F. W. Takes, W. A. Kosters,
82 "Computing the Eccentricity Distribution of Large Graphs."
83 Algorithms, 2013
84 https://www.mdpi.com/1999-4893/6/1/100
85 .. [3] M. Borassi, P. Crescenzi, M. Habib, W. A. Kosters, A. Marino, F. W. Takes,
86 "Fast diameter and radius BFS-based computation in (weakly connected) real-world graphs: With an application to the six degrees of separation games. "
87 Theoretical Computer Science, 2015
88 https://www.sciencedirect.com/science/article/pii/S0304397515001644
89 """
90 # init variables
91 degrees = dict(G.degree()) # start with the highest degree node
92 minlowernode = max(degrees, key=degrees.get)
93 N = len(degrees) # number of nodes
94 # alternate between smallest lower and largest upper bound
95 high = False
96 # status variables
97 ecc_lower = dict.fromkeys(G, 0)
98 ecc_upper = dict.fromkeys(G, N)
99 candidates = set(G)
101 # (re)set bound extremes
102 minlower = N
103 maxlower = 0
104 minupper = N
105 maxupper = 0
107 # repeat the following until there are no more candidates
108 while candidates:
109 if high:
110 current = maxuppernode # select node with largest upper bound
111 else:
112 current = minlowernode # select node with smallest lower bound
113 high = not high
115 # get distances from/to current node and derive eccentricity
116 dist = nx.shortest_path_length(G, source=current, weight=weight)
118 if len(dist) != N:
119 msg = "Cannot compute metric because graph is not connected."
120 raise nx.NetworkXError(msg)
121 current_ecc = max(dist.values())
123 # print status update
124 # print ("ecc of " + str(current) + " (" + str(ecc_lower[current]) + "/"
125 # + str(ecc_upper[current]) + ", deg: " + str(dist[current]) + ") is "
126 # + str(current_ecc))
127 # print(ecc_upper)
129 # (re)set bound extremes
130 maxuppernode = None
131 minlowernode = None
133 # update node bounds
134 for i in candidates:
135 # update eccentricity bounds
136 d = dist[i]
137 ecc_lower[i] = low = max(ecc_lower[i], max(d, (current_ecc - d)))
138 ecc_upper[i] = upp = min(ecc_upper[i], current_ecc + d)
140 # update min/max values of lower and upper bounds
141 minlower = min(ecc_lower[i], minlower)
142 maxlower = max(ecc_lower[i], maxlower)
143 minupper = min(ecc_upper[i], minupper)
144 maxupper = max(ecc_upper[i], maxupper)
146 # update candidate set
147 if compute == "diameter":
148 ruled_out = {
149 i
150 for i in candidates
151 if ecc_upper[i] <= maxlower and 2 * ecc_lower[i] >= maxupper
152 }
153 elif compute == "radius":
154 ruled_out = {
155 i
156 for i in candidates
157 if ecc_lower[i] >= minupper and ecc_upper[i] + 1 <= 2 * minlower
158 }
159 elif compute == "periphery":
160 ruled_out = {
161 i
162 for i in candidates
163 if ecc_upper[i] < maxlower
164 and (maxlower == maxupper or ecc_lower[i] > maxupper)
165 }
166 elif compute == "center":
167 ruled_out = {
168 i
169 for i in candidates
170 if ecc_lower[i] > minupper
171 and (minlower == minupper or ecc_upper[i] + 1 < 2 * minlower)
172 }
173 elif compute == "eccentricities":
174 ruled_out = set()
175 else:
176 msg = "compute must be one of 'diameter', 'radius', 'periphery', 'center', 'eccentricities'"
177 raise ValueError(msg)
179 ruled_out.update(i for i in candidates if ecc_lower[i] == ecc_upper[i])
180 candidates -= ruled_out
182 # for i in ruled_out:
183 # print("removing %g: ecc_u: %g maxl: %g ecc_l: %g maxu: %g"%
184 # (i,ecc_upper[i],maxlower,ecc_lower[i],maxupper))
185 # print("node %g: ecc_u: %g maxl: %g ecc_l: %g maxu: %g"%
186 # (4,ecc_upper[4],maxlower,ecc_lower[4],maxupper))
187 # print("NODE 4: %g"%(ecc_upper[4] <= maxlower))
188 # print("NODE 4: %g"%(2 * ecc_lower[4] >= maxupper))
189 # print("NODE 4: %g"%(ecc_upper[4] <= maxlower
190 # and 2 * ecc_lower[4] >= maxupper))
192 # updating maxuppernode and minlowernode for selection in next round
193 for i in candidates:
194 if (
195 minlowernode is None
196 or (
197 ecc_lower[i] == ecc_lower[minlowernode]
198 and degrees[i] > degrees[minlowernode]
199 )
200 or (ecc_lower[i] < ecc_lower[minlowernode])
201 ):
202 minlowernode = i
204 if (
205 maxuppernode is None
206 or (
207 ecc_upper[i] == ecc_upper[maxuppernode]
208 and degrees[i] > degrees[maxuppernode]
209 )
210 or (ecc_upper[i] > ecc_upper[maxuppernode])
211 ):
212 maxuppernode = i
214 # print status update
215 # print (" min=" + str(minlower) + "/" + str(minupper) +
216 # " max=" + str(maxlower) + "/" + str(maxupper) +
217 # " candidates: " + str(len(candidates)))
218 # print("cand:",candidates)
219 # print("ecc_l",ecc_lower)
220 # print("ecc_u",ecc_upper)
221 # wait = input("press Enter to continue")
223 # return the correct value of the requested metric
224 if compute == "diameter":
225 return maxlower
226 if compute == "radius":
227 return minupper
228 if compute == "periphery":
229 p = [v for v in G if ecc_lower[v] == maxlower]
230 return p
231 if compute == "center":
232 c = [v for v in G if ecc_upper[v] == minupper]
233 return c
234 if compute == "eccentricities":
235 return ecc_lower
236 return None
239@nx._dispatch(edge_attrs="weight")
240def eccentricity(G, v=None, sp=None, weight=None):
241 """Returns the eccentricity of nodes in G.
243 The eccentricity of a node v is the maximum distance from v to
244 all other nodes in G.
246 Parameters
247 ----------
248 G : NetworkX graph
249 A graph
251 v : node, optional
252 Return value of specified node
254 sp : dict of dicts, optional
255 All pairs shortest path lengths as a dictionary of dictionaries
257 weight : string, function, or None (default=None)
258 If this is a string, then edge weights will be accessed via the
259 edge attribute with this key (that is, the weight of the edge
260 joining `u` to `v` will be ``G.edges[u, v][weight]``). If no
261 such edge attribute exists, the weight of the edge is assumed to
262 be one.
264 If this is a function, the weight of an edge is the value
265 returned by the function. The function must accept exactly three
266 positional arguments: the two endpoints of an edge and the
267 dictionary of edge attributes for that edge. The function must
268 return a number.
270 If this is None, every edge has weight/distance/cost 1.
272 Weights stored as floating point values can lead to small round-off
273 errors in distances. Use integer weights to avoid this.
275 Weights should be positive, since they are distances.
277 Returns
278 -------
279 ecc : dictionary
280 A dictionary of eccentricity values keyed by node.
282 Examples
283 --------
284 >>> G = nx.Graph([(1, 2), (1, 3), (1, 4), (3, 4), (3, 5), (4, 5)])
285 >>> dict(nx.eccentricity(G))
286 {1: 2, 2: 3, 3: 2, 4: 2, 5: 3}
288 >>> dict(nx.eccentricity(G, v=[1, 5])) # This returns the eccentricity of node 1 & 5
289 {1: 2, 5: 3}
291 """
292 # if v is None: # none, use entire graph
293 # nodes=G.nodes()
294 # elif v in G: # is v a single node
295 # nodes=[v]
296 # else: # assume v is a container of nodes
297 # nodes=v
298 order = G.order()
299 e = {}
300 for n in G.nbunch_iter(v):
301 if sp is None:
302 length = nx.shortest_path_length(G, source=n, weight=weight)
304 L = len(length)
305 else:
306 try:
307 length = sp[n]
308 L = len(length)
309 except TypeError as err:
310 raise nx.NetworkXError('Format of "sp" is invalid.') from err
311 if L != order:
312 if G.is_directed():
313 msg = (
314 "Found infinite path length because the digraph is not"
315 " strongly connected"
316 )
317 else:
318 msg = "Found infinite path length because the graph is not" " connected"
319 raise nx.NetworkXError(msg)
321 e[n] = max(length.values())
323 if v in G:
324 return e[v] # return single value
325 return e
328@nx._dispatch(edge_attrs="weight")
329def diameter(G, e=None, usebounds=False, weight=None):
330 """Returns the diameter of the graph G.
332 The diameter is the maximum eccentricity.
334 Parameters
335 ----------
336 G : NetworkX graph
337 A graph
339 e : eccentricity dictionary, optional
340 A precomputed dictionary of eccentricities.
342 weight : string, function, or None
343 If this is a string, then edge weights will be accessed via the
344 edge attribute with this key (that is, the weight of the edge
345 joining `u` to `v` will be ``G.edges[u, v][weight]``). If no
346 such edge attribute exists, the weight of the edge is assumed to
347 be one.
349 If this is a function, the weight of an edge is the value
350 returned by the function. The function must accept exactly three
351 positional arguments: the two endpoints of an edge and the
352 dictionary of edge attributes for that edge. The function must
353 return a number.
355 If this is None, every edge has weight/distance/cost 1.
357 Weights stored as floating point values can lead to small round-off
358 errors in distances. Use integer weights to avoid this.
360 Weights should be positive, since they are distances.
362 Returns
363 -------
364 d : integer
365 Diameter of graph
367 Examples
368 --------
369 >>> G = nx.Graph([(1, 2), (1, 3), (1, 4), (3, 4), (3, 5), (4, 5)])
370 >>> nx.diameter(G)
371 3
373 See Also
374 --------
375 eccentricity
376 """
377 if usebounds is True and e is None and not G.is_directed():
378 return _extrema_bounding(G, compute="diameter", weight=weight)
379 if e is None:
380 e = eccentricity(G, weight=weight)
381 return max(e.values())
384@nx._dispatch(edge_attrs="weight")
385def periphery(G, e=None, usebounds=False, weight=None):
386 """Returns the periphery of the graph G.
388 The periphery is the set of nodes with eccentricity equal to the diameter.
390 Parameters
391 ----------
392 G : NetworkX graph
393 A graph
395 e : eccentricity dictionary, optional
396 A precomputed dictionary of eccentricities.
398 weight : string, function, or None
399 If this is a string, then edge weights will be accessed via the
400 edge attribute with this key (that is, the weight of the edge
401 joining `u` to `v` will be ``G.edges[u, v][weight]``). If no
402 such edge attribute exists, the weight of the edge is assumed to
403 be one.
405 If this is a function, the weight of an edge is the value
406 returned by the function. The function must accept exactly three
407 positional arguments: the two endpoints of an edge and the
408 dictionary of edge attributes for that edge. The function must
409 return a number.
411 If this is None, every edge has weight/distance/cost 1.
413 Weights stored as floating point values can lead to small round-off
414 errors in distances. Use integer weights to avoid this.
416 Weights should be positive, since they are distances.
418 Returns
419 -------
420 p : list
421 List of nodes in periphery
423 Examples
424 --------
425 >>> G = nx.Graph([(1, 2), (1, 3), (1, 4), (3, 4), (3, 5), (4, 5)])
426 >>> nx.periphery(G)
427 [2, 5]
429 See Also
430 --------
431 barycenter
432 center
433 """
434 if usebounds is True and e is None and not G.is_directed():
435 return _extrema_bounding(G, compute="periphery", weight=weight)
436 if e is None:
437 e = eccentricity(G, weight=weight)
438 diameter = max(e.values())
439 p = [v for v in e if e[v] == diameter]
440 return p
443@nx._dispatch(edge_attrs="weight")
444def radius(G, e=None, usebounds=False, weight=None):
445 """Returns the radius of the graph G.
447 The radius is the minimum eccentricity.
449 Parameters
450 ----------
451 G : NetworkX graph
452 A graph
454 e : eccentricity dictionary, optional
455 A precomputed dictionary of eccentricities.
457 weight : string, function, or None
458 If this is a string, then edge weights will be accessed via the
459 edge attribute with this key (that is, the weight of the edge
460 joining `u` to `v` will be ``G.edges[u, v][weight]``). If no
461 such edge attribute exists, the weight of the edge is assumed to
462 be one.
464 If this is a function, the weight of an edge is the value
465 returned by the function. The function must accept exactly three
466 positional arguments: the two endpoints of an edge and the
467 dictionary of edge attributes for that edge. The function must
468 return a number.
470 If this is None, every edge has weight/distance/cost 1.
472 Weights stored as floating point values can lead to small round-off
473 errors in distances. Use integer weights to avoid this.
475 Weights should be positive, since they are distances.
477 Returns
478 -------
479 r : integer
480 Radius of graph
482 Examples
483 --------
484 >>> G = nx.Graph([(1, 2), (1, 3), (1, 4), (3, 4), (3, 5), (4, 5)])
485 >>> nx.radius(G)
486 2
488 """
489 if usebounds is True and e is None and not G.is_directed():
490 return _extrema_bounding(G, compute="radius", weight=weight)
491 if e is None:
492 e = eccentricity(G, weight=weight)
493 return min(e.values())
496@nx._dispatch(edge_attrs="weight")
497def center(G, e=None, usebounds=False, weight=None):
498 """Returns the center of the graph G.
500 The center is the set of nodes with eccentricity equal to radius.
502 Parameters
503 ----------
504 G : NetworkX graph
505 A graph
507 e : eccentricity dictionary, optional
508 A precomputed dictionary of eccentricities.
510 weight : string, function, or None
511 If this is a string, then edge weights will be accessed via the
512 edge attribute with this key (that is, the weight of the edge
513 joining `u` to `v` will be ``G.edges[u, v][weight]``). If no
514 such edge attribute exists, the weight of the edge is assumed to
515 be one.
517 If this is a function, the weight of an edge is the value
518 returned by the function. The function must accept exactly three
519 positional arguments: the two endpoints of an edge and the
520 dictionary of edge attributes for that edge. The function must
521 return a number.
523 If this is None, every edge has weight/distance/cost 1.
525 Weights stored as floating point values can lead to small round-off
526 errors in distances. Use integer weights to avoid this.
528 Weights should be positive, since they are distances.
530 Returns
531 -------
532 c : list
533 List of nodes in center
535 Examples
536 --------
537 >>> G = nx.Graph([(1, 2), (1, 3), (1, 4), (3, 4), (3, 5), (4, 5)])
538 >>> list(nx.center(G))
539 [1, 3, 4]
541 See Also
542 --------
543 barycenter
544 periphery
545 """
546 if usebounds is True and e is None and not G.is_directed():
547 return _extrema_bounding(G, compute="center", weight=weight)
548 if e is None:
549 e = eccentricity(G, weight=weight)
550 radius = min(e.values())
551 p = [v for v in e if e[v] == radius]
552 return p
555@nx._dispatch(edge_attrs="weight")
556def barycenter(G, weight=None, attr=None, sp=None):
557 r"""Calculate barycenter of a connected graph, optionally with edge weights.
559 The :dfn:`barycenter` a
560 :func:`connected <networkx.algorithms.components.is_connected>` graph
561 :math:`G` is the subgraph induced by the set of its nodes :math:`v`
562 minimizing the objective function
564 .. math::
566 \sum_{u \in V(G)} d_G(u, v),
568 where :math:`d_G` is the (possibly weighted) :func:`path length
569 <networkx.algorithms.shortest_paths.generic.shortest_path_length>`.
570 The barycenter is also called the :dfn:`median`. See [West01]_, p. 78.
572 Parameters
573 ----------
574 G : :class:`networkx.Graph`
575 The connected graph :math:`G`.
576 weight : :class:`str`, optional
577 Passed through to
578 :func:`~networkx.algorithms.shortest_paths.generic.shortest_path_length`.
579 attr : :class:`str`, optional
580 If given, write the value of the objective function to each node's
581 `attr` attribute. Otherwise do not store the value.
582 sp : dict of dicts, optional
583 All pairs shortest path lengths as a dictionary of dictionaries
585 Returns
586 -------
587 list
588 Nodes of `G` that induce the barycenter of `G`.
590 Raises
591 ------
592 NetworkXNoPath
593 If `G` is disconnected. `G` may appear disconnected to
594 :func:`barycenter` if `sp` is given but is missing shortest path
595 lengths for any pairs.
596 ValueError
597 If `sp` and `weight` are both given.
599 Examples
600 --------
601 >>> G = nx.Graph([(1, 2), (1, 3), (1, 4), (3, 4), (3, 5), (4, 5)])
602 >>> nx.barycenter(G)
603 [1, 3, 4]
605 See Also
606 --------
607 center
608 periphery
609 """
610 if sp is None:
611 sp = nx.shortest_path_length(G, weight=weight)
612 else:
613 sp = sp.items()
614 if weight is not None:
615 raise ValueError("Cannot use both sp, weight arguments together")
616 smallest, barycenter_vertices, n = float("inf"), [], len(G)
617 for v, dists in sp:
618 if len(dists) < n:
619 raise nx.NetworkXNoPath(
620 f"Input graph {G} is disconnected, so every induced subgraph "
621 "has infinite barycentricity."
622 )
623 barycentricity = sum(dists.values())
624 if attr is not None:
625 G.nodes[v][attr] = barycentricity
626 if barycentricity < smallest:
627 smallest = barycentricity
628 barycenter_vertices = [v]
629 elif barycentricity == smallest:
630 barycenter_vertices.append(v)
631 return barycenter_vertices
634def _count_lu_permutations(perm_array):
635 """Counts the number of permutations in SuperLU perm_c or perm_r"""
636 perm_cnt = 0
637 arr = perm_array.tolist()
638 for i in range(len(arr)):
639 if i != arr[i]:
640 perm_cnt += 1
641 n = arr.index(i)
642 arr[n] = arr[i]
643 arr[i] = i
645 return perm_cnt
648@not_implemented_for("directed")
649@nx._dispatch(edge_attrs="weight")
650def resistance_distance(G, nodeA=None, nodeB=None, weight=None, invert_weight=True):
651 """Returns the resistance distance between every pair of nodes on graph G.
653 The resistance distance between two nodes of a graph is akin to treating
654 the graph as a grid of resistors with a resistance equal to the provided
655 weight [1]_, [2]_.
657 If weight is not provided, then a weight of 1 is used for all edges.
659 If two nodes are the same, the resistance distance is zero.
661 Parameters
662 ----------
663 G : NetworkX graph
664 A graph
666 nodeA : node or None, optional (default=None)
667 A node within graph G.
668 If None, compute resistance distance using all nodes as source nodes.
670 nodeB : node or None, optional (default=None)
671 A node within graph G.
672 If None, compute resistance distance using all nodes as target nodes.
674 weight : string or None, optional (default=None)
675 The edge data key used to compute the resistance distance.
676 If None, then each edge has weight 1.
678 invert_weight : boolean (default=True)
679 Proper calculation of resistance distance requires building the
680 Laplacian matrix with the reciprocal of the weight. Not required
681 if the weight is already inverted. Weight cannot be zero.
683 Returns
684 -------
685 rd : dict or float
686 If `nodeA` and `nodeB` are given, resistance distance between `nodeA`
687 and `nodeB`. If `nodeA` or `nodeB` is unspecified (the default), a
688 dictionary of nodes with resistance distances as the value.
690 Raises
691 ------
692 NetworkXNotImplemented
693 If `G` is a directed graph.
695 NetworkXError
696 If `G` is not connected, or contains no nodes,
697 or `nodeA` is not in `G` or `nodeB` is not in `G`.
699 Examples
700 --------
701 >>> G = nx.Graph([(1, 2), (1, 3), (1, 4), (3, 4), (3, 5), (4, 5)])
702 >>> round(nx.resistance_distance(G, 1, 3), 10)
703 0.625
705 Notes
706 -----
707 The implementation is based on Theorem A in [2]_. Self-loops are ignored.
708 Multi-edges are contracted in one edge with weight equal to the harmonic sum of the weights.
710 References
711 ----------
712 .. [1] Wikipedia
713 "Resistance distance."
714 https://en.wikipedia.org/wiki/Resistance_distance
715 .. [2] D. J. Klein and M. Randic.
716 Resistance distance.
717 J. of Math. Chem. 12:81-95, 1993.
718 """
719 import numpy as np
721 if len(G) == 0:
722 raise nx.NetworkXError("Graph G must contain at least one node.")
723 if not nx.is_connected(G):
724 raise nx.NetworkXError("Graph G must be strongly connected.")
725 if nodeA is not None and nodeA not in G:
726 raise nx.NetworkXError("Node A is not in graph G.")
727 if nodeB is not None and nodeB not in G:
728 raise nx.NetworkXError("Node B is not in graph G.")
730 G = G.copy()
731 node_list = list(G)
733 # Invert weights
734 if invert_weight and weight is not None:
735 if G.is_multigraph():
736 for u, v, k, d in G.edges(keys=True, data=True):
737 d[weight] = 1 / d[weight]
738 else:
739 for u, v, d in G.edges(data=True):
740 d[weight] = 1 / d[weight]
742 # Compute resistance distance using the Pseudo-inverse of the Laplacian
743 # Self-loops are ignored
744 L = nx.laplacian_matrix(G, weight=weight).todense()
745 Linv = np.linalg.pinv(L, hermitian=True)
747 # Return relevant distances
748 if nodeA is not None and nodeB is not None:
749 i = node_list.index(nodeA)
750 j = node_list.index(nodeB)
751 return Linv[i, i] + Linv[j, j] - Linv[i, j] - Linv[j, i]
753 elif nodeA is not None:
754 i = node_list.index(nodeA)
755 d = {}
756 for n in G:
757 j = node_list.index(n)
758 d[n] = Linv[i, i] + Linv[j, j] - Linv[i, j] - Linv[j, i]
759 return d
761 elif nodeB is not None:
762 j = node_list.index(nodeB)
763 d = {}
764 for n in G:
765 i = node_list.index(n)
766 d[n] = Linv[i, i] + Linv[j, j] - Linv[i, j] - Linv[j, i]
767 return d
769 else:
770 d = {}
771 for n in G:
772 i = node_list.index(n)
773 d[n] = {}
774 for n2 in G:
775 j = node_list.index(n2)
776 d[n][n2] = Linv[i, i] + Linv[j, j] - Linv[i, j] - Linv[j, i]
777 return d
778 # Replace with collapsing topology or approximated zero?
780 # Using determinants to compute the effective resistance is more memory
781 # efficient than directly calculating the pseudo-inverse
782 L = nx.laplacian_matrix(G, node_list, weight=weight).asformat("csc")
783 indices = list(range(L.shape[0]))
784 # w/ nodeA removed
785 indices.remove(node_list.index(nodeA))
786 L_a = L[indices, :][:, indices]
787 # Both nodeA and nodeB removed
788 indices.remove(node_list.index(nodeB))
789 L_ab = L[indices, :][:, indices]
791 # Factorize Laplacian submatrixes and extract diagonals
792 # Order the diagonals to minimize the likelihood over overflows
793 # during computing the determinant
794 lu_a = sp.sparse.linalg.splu(L_a, options={"SymmetricMode": True})
795 LdiagA = lu_a.U.diagonal()
796 LdiagA_s = np.prod(np.sign(LdiagA)) * np.prod(lu_a.L.diagonal())
797 LdiagA_s *= (-1) ** _count_lu_permutations(lu_a.perm_r)
798 LdiagA_s *= (-1) ** _count_lu_permutations(lu_a.perm_c)
799 LdiagA = np.absolute(LdiagA)
800 LdiagA = np.sort(LdiagA)
802 lu_ab = sp.sparse.linalg.splu(L_ab, options={"SymmetricMode": True})
803 LdiagAB = lu_ab.U.diagonal()
804 LdiagAB_s = np.prod(np.sign(LdiagAB)) * np.prod(lu_ab.L.diagonal())
805 LdiagAB_s *= (-1) ** _count_lu_permutations(lu_ab.perm_r)
806 LdiagAB_s *= (-1) ** _count_lu_permutations(lu_ab.perm_c)
807 LdiagAB = np.absolute(LdiagAB)
808 LdiagAB = np.sort(LdiagAB)
810 # Calculate the ratio of determinant, rd = det(L_ab)/det(L_a)
811 Ldet = np.prod(np.divide(np.append(LdiagAB, [1]), LdiagA))
812 rd = Ldet * LdiagAB_s / LdiagA_s
814 return rd
817@nx.utils.not_implemented_for("directed")
818@nx._dispatch(edge_attrs="weight")
819def kemeny_constant(G, *, weight=None):
820 """Returns the Kemeny constant of the given graph.
822 The *Kemeny constant* (or Kemeny's constant) of a graph `G`
823 can be computed by regarding the graph as a Markov chain.
824 The Kemeny constant is then the expected number of time steps
825 to transition from a starting state i to a random destination state
826 sampled from the Markov chain's stationary distribution.
827 The Kemeny constant is independent of the chosen initial state [1]_.
829 The Kemeny constant measures the time needed for spreading
830 across a graph. Low values indicate a closely connected graph
831 whereas high values indicate a spread-out graph.
833 If weight is not provided, then a weight of 1 is used for all edges.
835 Since `G` represents a Markov chain, the weights must be positive.
837 Parameters
838 ----------
839 G : NetworkX graph
841 weight : string or None, optional (default=None)
842 The edge data key used to compute the Kemeny constant.
843 If None, then each edge has weight 1.
845 Returns
846 -------
847 K : float
848 The Kemeny constant of the graph `G`.
850 Raises
851 ------
852 NetworkXNotImplemented
853 If the graph `G` is directed.
855 NetworkXError
856 If the graph `G` is not connected, or contains no nodes,
857 or has edges with negative weights.
859 Examples
860 --------
861 >>> G = nx.complete_graph(5)
862 >>> round(nx.kemeny_constant(G), 10)
863 3.2
865 Notes
866 -----
867 The implementation is based on equation (3.3) in [2]_.
868 Self-loops are allowed and indicate a Markov chain where
869 the state can remain the same. Multi-edges are contracted
870 in one edge with weight equal to the sum of the weights.
872 References
873 ----------
874 .. [1] Wikipedia
875 "Kemeny's constant."
876 https://en.wikipedia.org/wiki/Kemeny%27s_constant
877 .. [2] Lovász L.
878 Random walks on graphs: A survey.
879 Paul Erdös is Eighty, vol. 2, Bolyai Society,
880 Mathematical Studies, Keszthely, Hungary (1993), pp. 1-46
881 """
882 import numpy as np
883 import scipy as sp
885 if len(G) == 0:
886 raise nx.NetworkXError("Graph G must contain at least one node.")
887 if not nx.is_connected(G):
888 raise nx.NetworkXError("Graph G must be connected.")
889 if nx.is_negatively_weighted(G, weight=weight):
890 raise nx.NetworkXError("The weights of graph G must be nonnegative.")
892 # Compute matrix H = D^-1/2 A D^-1/2
893 A = nx.adjacency_matrix(G, weight=weight)
894 n, m = A.shape
895 diags = A.sum(axis=1)
896 with np.errstate(divide="ignore"):
897 diags_sqrt = 1.0 / np.sqrt(diags)
898 diags_sqrt[np.isinf(diags_sqrt)] = 0
899 DH = sp.sparse.csr_array(sp.sparse.spdiags(diags_sqrt, 0, m, n, format="csr"))
900 H = DH @ (A @ DH)
902 # Compute eigenvalues of H
903 eig = np.sort(sp.linalg.eigvalsh(H.todense()))
905 # Compute the Kemeny constant
906 return np.sum(1 / (1 - eig[:-1]))