matching.py 21 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585
  1. # This module uses material from the Wikipedia article Hopcroft--Karp algorithm
  2. # <https://en.wikipedia.org/wiki/Hopcroft%E2%80%93Karp_algorithm>, accessed on
  3. # January 3, 2015, which is released under the Creative Commons
  4. # Attribution-Share-Alike License 3.0
  5. # <http://creativecommons.org/licenses/by-sa/3.0/>. That article includes
  6. # pseudocode, which has been translated into the corresponding Python code.
  7. #
  8. # Portions of this module use code from David Eppstein's Python Algorithms and
  9. # Data Structures (PADS) library, which is dedicated to the public domain (for
  10. # proof, see <http://www.ics.uci.edu/~eppstein/PADS/ABOUT-PADS.txt>).
  11. """Provides functions for computing maximum cardinality matchings and minimum
  12. weight full matchings in a bipartite graph.
  13. If you don't care about the particular implementation of the maximum matching
  14. algorithm, simply use the :func:`maximum_matching`. If you do care, you can
  15. import one of the named maximum matching algorithms directly.
  16. For example, to find a maximum matching in the complete bipartite graph with
  17. two vertices on the left and three vertices on the right:
  18. >>> G = nx.complete_bipartite_graph(2, 3)
  19. >>> left, right = nx.bipartite.sets(G)
  20. >>> list(left)
  21. [0, 1]
  22. >>> list(right)
  23. [2, 3, 4]
  24. >>> nx.bipartite.maximum_matching(G)
  25. {0: 2, 1: 3, 2: 0, 3: 1}
  26. The dictionary returned by :func:`maximum_matching` includes a mapping for
  27. vertices in both the left and right vertex sets.
  28. Similarly, :func:`minimum_weight_full_matching` produces, for a complete
  29. weighted bipartite graph, a matching whose cardinality is the cardinality of
  30. the smaller of the two partitions, and for which the sum of the weights of the
  31. edges included in the matching is minimal.
  32. """
  33. import collections
  34. import itertools
  35. import networkx as nx
  36. from networkx.algorithms.bipartite import sets as bipartite_sets
  37. from networkx.algorithms.bipartite.matrix import biadjacency_matrix
  38. __all__ = [
  39. "maximum_matching",
  40. "hopcroft_karp_matching",
  41. "eppstein_matching",
  42. "to_vertex_cover",
  43. "minimum_weight_full_matching",
  44. ]
  45. INFINITY = float("inf")
  46. def hopcroft_karp_matching(G, top_nodes=None):
  47. """Returns the maximum cardinality matching of the bipartite graph `G`.
  48. A matching is a set of edges that do not share any nodes. A maximum
  49. cardinality matching is a matching with the most edges possible. It
  50. is not always unique. Finding a matching in a bipartite graph can be
  51. treated as a networkx flow problem.
  52. The functions ``hopcroft_karp_matching`` and ``maximum_matching``
  53. are aliases of the same function.
  54. Parameters
  55. ----------
  56. G : NetworkX graph
  57. Undirected bipartite graph
  58. top_nodes : container of nodes
  59. Container with all nodes in one bipartite node set. If not supplied
  60. it will be computed. But if more than one solution exists an exception
  61. will be raised.
  62. Returns
  63. -------
  64. matches : dictionary
  65. The matching is returned as a dictionary, `matches`, such that
  66. ``matches[v] == w`` if node `v` is matched to node `w`. Unmatched
  67. nodes do not occur as a key in `matches`.
  68. Raises
  69. ------
  70. AmbiguousSolution
  71. Raised if the input bipartite graph is disconnected and no container
  72. with all nodes in one bipartite set is provided. When determining
  73. the nodes in each bipartite set more than one valid solution is
  74. possible if the input graph is disconnected.
  75. Notes
  76. -----
  77. This function is implemented with the `Hopcroft--Karp matching algorithm
  78. <https://en.wikipedia.org/wiki/Hopcroft%E2%80%93Karp_algorithm>`_ for
  79. bipartite graphs.
  80. See :mod:`bipartite documentation <networkx.algorithms.bipartite>`
  81. for further details on how bipartite graphs are handled in NetworkX.
  82. See Also
  83. --------
  84. maximum_matching
  85. hopcroft_karp_matching
  86. eppstein_matching
  87. References
  88. ----------
  89. .. [1] John E. Hopcroft and Richard M. Karp. "An n^{5 / 2} Algorithm for
  90. Maximum Matchings in Bipartite Graphs" In: **SIAM Journal of Computing**
  91. 2.4 (1973), pp. 225--231. <https://doi.org/10.1137/0202019>.
  92. """
  93. # First we define some auxiliary search functions.
  94. #
  95. # If you are a human reading these auxiliary search functions, the "global"
  96. # variables `leftmatches`, `rightmatches`, `distances`, etc. are defined
  97. # below the functions, so that they are initialized close to the initial
  98. # invocation of the search functions.
  99. def breadth_first_search():
  100. for v in left:
  101. if leftmatches[v] is None:
  102. distances[v] = 0
  103. queue.append(v)
  104. else:
  105. distances[v] = INFINITY
  106. distances[None] = INFINITY
  107. while queue:
  108. v = queue.popleft()
  109. if distances[v] < distances[None]:
  110. for u in G[v]:
  111. if distances[rightmatches[u]] is INFINITY:
  112. distances[rightmatches[u]] = distances[v] + 1
  113. queue.append(rightmatches[u])
  114. return distances[None] is not INFINITY
  115. def depth_first_search(v):
  116. if v is not None:
  117. for u in G[v]:
  118. if distances[rightmatches[u]] == distances[v] + 1:
  119. if depth_first_search(rightmatches[u]):
  120. rightmatches[u] = v
  121. leftmatches[v] = u
  122. return True
  123. distances[v] = INFINITY
  124. return False
  125. return True
  126. # Initialize the "global" variables that maintain state during the search.
  127. left, right = bipartite_sets(G, top_nodes)
  128. leftmatches = {v: None for v in left}
  129. rightmatches = {v: None for v in right}
  130. distances = {}
  131. queue = collections.deque()
  132. # Implementation note: this counter is incremented as pairs are matched but
  133. # it is currently not used elsewhere in the computation.
  134. num_matched_pairs = 0
  135. while breadth_first_search():
  136. for v in left:
  137. if leftmatches[v] is None:
  138. if depth_first_search(v):
  139. num_matched_pairs += 1
  140. # Strip the entries matched to `None`.
  141. leftmatches = {k: v for k, v in leftmatches.items() if v is not None}
  142. rightmatches = {k: v for k, v in rightmatches.items() if v is not None}
  143. # At this point, the left matches and the right matches are inverses of one
  144. # another. In other words,
  145. #
  146. # leftmatches == {v, k for k, v in rightmatches.items()}
  147. #
  148. # Finally, we combine both the left matches and right matches.
  149. return dict(itertools.chain(leftmatches.items(), rightmatches.items()))
  150. def eppstein_matching(G, top_nodes=None):
  151. """Returns the maximum cardinality matching of the bipartite graph `G`.
  152. Parameters
  153. ----------
  154. G : NetworkX graph
  155. Undirected bipartite graph
  156. top_nodes : container
  157. Container with all nodes in one bipartite node set. If not supplied
  158. it will be computed. But if more than one solution exists an exception
  159. will be raised.
  160. Returns
  161. -------
  162. matches : dictionary
  163. The matching is returned as a dictionary, `matching`, such that
  164. ``matching[v] == w`` if node `v` is matched to node `w`. Unmatched
  165. nodes do not occur as a key in `matching`.
  166. Raises
  167. ------
  168. AmbiguousSolution
  169. Raised if the input bipartite graph is disconnected and no container
  170. with all nodes in one bipartite set is provided. When determining
  171. the nodes in each bipartite set more than one valid solution is
  172. possible if the input graph is disconnected.
  173. Notes
  174. -----
  175. This function is implemented with David Eppstein's version of the algorithm
  176. Hopcroft--Karp algorithm (see :func:`hopcroft_karp_matching`), which
  177. originally appeared in the `Python Algorithms and Data Structures library
  178. (PADS) <http://www.ics.uci.edu/~eppstein/PADS/ABOUT-PADS.txt>`_.
  179. See :mod:`bipartite documentation <networkx.algorithms.bipartite>`
  180. for further details on how bipartite graphs are handled in NetworkX.
  181. See Also
  182. --------
  183. hopcroft_karp_matching
  184. """
  185. # Due to its original implementation, a directed graph is needed
  186. # so that the two sets of bipartite nodes can be distinguished
  187. left, right = bipartite_sets(G, top_nodes)
  188. G = nx.DiGraph(G.edges(left))
  189. # initialize greedy matching (redundant, but faster than full search)
  190. matching = {}
  191. for u in G:
  192. for v in G[u]:
  193. if v not in matching:
  194. matching[v] = u
  195. break
  196. while True:
  197. # structure residual graph into layers
  198. # pred[u] gives the neighbor in the previous layer for u in U
  199. # preds[v] gives a list of neighbors in the previous layer for v in V
  200. # unmatched gives a list of unmatched vertices in final layer of V,
  201. # and is also used as a flag value for pred[u] when u is in the first
  202. # layer
  203. preds = {}
  204. unmatched = []
  205. pred = {u: unmatched for u in G}
  206. for v in matching:
  207. del pred[matching[v]]
  208. layer = list(pred)
  209. # repeatedly extend layering structure by another pair of layers
  210. while layer and not unmatched:
  211. newLayer = {}
  212. for u in layer:
  213. for v in G[u]:
  214. if v not in preds:
  215. newLayer.setdefault(v, []).append(u)
  216. layer = []
  217. for v in newLayer:
  218. preds[v] = newLayer[v]
  219. if v in matching:
  220. layer.append(matching[v])
  221. pred[matching[v]] = v
  222. else:
  223. unmatched.append(v)
  224. # did we finish layering without finding any alternating paths?
  225. if not unmatched:
  226. # TODO - The lines between --- were unused and were thus commented
  227. # out. This whole commented chunk should be reviewed to determine
  228. # whether it should be built upon or completely removed.
  229. # ---
  230. # unlayered = {}
  231. # for u in G:
  232. # # TODO Why is extra inner loop necessary?
  233. # for v in G[u]:
  234. # if v not in preds:
  235. # unlayered[v] = None
  236. # ---
  237. # TODO Originally, this function returned a three-tuple:
  238. #
  239. # return (matching, list(pred), list(unlayered))
  240. #
  241. # For some reason, the documentation for this function
  242. # indicated that the second and third elements of the returned
  243. # three-tuple would be the vertices in the left and right vertex
  244. # sets, respectively, that are also in the maximum independent set.
  245. # However, what I think the author meant was that the second
  246. # element is the list of vertices that were unmatched and the third
  247. # element was the list of vertices that were matched. Since that
  248. # seems to be the case, they don't really need to be returned,
  249. # since that information can be inferred from the matching
  250. # dictionary.
  251. # All the matched nodes must be a key in the dictionary
  252. for key in matching.copy():
  253. matching[matching[key]] = key
  254. return matching
  255. # recursively search backward through layers to find alternating paths
  256. # recursion returns true if found path, false otherwise
  257. def recurse(v):
  258. if v in preds:
  259. L = preds.pop(v)
  260. for u in L:
  261. if u in pred:
  262. pu = pred.pop(u)
  263. if pu is unmatched or recurse(pu):
  264. matching[v] = u
  265. return True
  266. return False
  267. for v in unmatched:
  268. recurse(v)
  269. def _is_connected_by_alternating_path(G, v, matched_edges, unmatched_edges, targets):
  270. """Returns True if and only if the vertex `v` is connected to one of
  271. the target vertices by an alternating path in `G`.
  272. An *alternating path* is a path in which every other edge is in the
  273. specified maximum matching (and the remaining edges in the path are not in
  274. the matching). An alternating path may have matched edges in the even
  275. positions or in the odd positions, as long as the edges alternate between
  276. 'matched' and 'unmatched'.
  277. `G` is an undirected bipartite NetworkX graph.
  278. `v` is a vertex in `G`.
  279. `matched_edges` is a set of edges present in a maximum matching in `G`.
  280. `unmatched_edges` is a set of edges not present in a maximum
  281. matching in `G`.
  282. `targets` is a set of vertices.
  283. """
  284. def _alternating_dfs(u, along_matched=True):
  285. """Returns True if and only if `u` is connected to one of the
  286. targets by an alternating path.
  287. `u` is a vertex in the graph `G`.
  288. If `along_matched` is True, this step of the depth-first search
  289. will continue only through edges in the given matching. Otherwise, it
  290. will continue only through edges *not* in the given matching.
  291. """
  292. visited = set()
  293. # Follow matched edges when depth is even,
  294. # and follow unmatched edges when depth is odd.
  295. initial_depth = 0 if along_matched else 1
  296. stack = [(u, iter(G[u]), initial_depth)]
  297. while stack:
  298. parent, children, depth = stack[-1]
  299. valid_edges = matched_edges if depth % 2 else unmatched_edges
  300. try:
  301. child = next(children)
  302. if child not in visited:
  303. if (parent, child) in valid_edges or (child, parent) in valid_edges:
  304. if child in targets:
  305. return True
  306. visited.add(child)
  307. stack.append((child, iter(G[child]), depth + 1))
  308. except StopIteration:
  309. stack.pop()
  310. return False
  311. # Check for alternating paths starting with edges in the matching, then
  312. # check for alternating paths starting with edges not in the
  313. # matching.
  314. return _alternating_dfs(v, along_matched=True) or _alternating_dfs(
  315. v, along_matched=False
  316. )
  317. def _connected_by_alternating_paths(G, matching, targets):
  318. """Returns the set of vertices that are connected to one of the target
  319. vertices by an alternating path in `G` or are themselves a target.
  320. An *alternating path* is a path in which every other edge is in the
  321. specified maximum matching (and the remaining edges in the path are not in
  322. the matching). An alternating path may have matched edges in the even
  323. positions or in the odd positions, as long as the edges alternate between
  324. 'matched' and 'unmatched'.
  325. `G` is an undirected bipartite NetworkX graph.
  326. `matching` is a dictionary representing a maximum matching in `G`, as
  327. returned by, for example, :func:`maximum_matching`.
  328. `targets` is a set of vertices.
  329. """
  330. # Get the set of matched edges and the set of unmatched edges. Only include
  331. # one version of each undirected edge (for example, include edge (1, 2) but
  332. # not edge (2, 1)). Using frozensets as an intermediary step we do not
  333. # require nodes to be orderable.
  334. edge_sets = {frozenset((u, v)) for u, v in matching.items()}
  335. matched_edges = {tuple(edge) for edge in edge_sets}
  336. unmatched_edges = {
  337. (u, v) for (u, v) in G.edges() if frozenset((u, v)) not in edge_sets
  338. }
  339. return {
  340. v
  341. for v in G
  342. if v in targets
  343. or _is_connected_by_alternating_path(
  344. G, v, matched_edges, unmatched_edges, targets
  345. )
  346. }
  347. def to_vertex_cover(G, matching, top_nodes=None):
  348. """Returns the minimum vertex cover corresponding to the given maximum
  349. matching of the bipartite graph `G`.
  350. Parameters
  351. ----------
  352. G : NetworkX graph
  353. Undirected bipartite graph
  354. matching : dictionary
  355. A dictionary whose keys are vertices in `G` and whose values are the
  356. distinct neighbors comprising the maximum matching for `G`, as returned
  357. by, for example, :func:`maximum_matching`. The dictionary *must*
  358. represent the maximum matching.
  359. top_nodes : container
  360. Container with all nodes in one bipartite node set. If not supplied
  361. it will be computed. But if more than one solution exists an exception
  362. will be raised.
  363. Returns
  364. -------
  365. vertex_cover : :class:`set`
  366. The minimum vertex cover in `G`.
  367. Raises
  368. ------
  369. AmbiguousSolution
  370. Raised if the input bipartite graph is disconnected and no container
  371. with all nodes in one bipartite set is provided. When determining
  372. the nodes in each bipartite set more than one valid solution is
  373. possible if the input graph is disconnected.
  374. Notes
  375. -----
  376. This function is implemented using the procedure guaranteed by `Konig's
  377. theorem
  378. <https://en.wikipedia.org/wiki/K%C3%B6nig%27s_theorem_%28graph_theory%29>`_,
  379. which proves an equivalence between a maximum matching and a minimum vertex
  380. cover in bipartite graphs.
  381. Since a minimum vertex cover is the complement of a maximum independent set
  382. for any graph, one can compute the maximum independent set of a bipartite
  383. graph this way:
  384. >>> G = nx.complete_bipartite_graph(2, 3)
  385. >>> matching = nx.bipartite.maximum_matching(G)
  386. >>> vertex_cover = nx.bipartite.to_vertex_cover(G, matching)
  387. >>> independent_set = set(G) - vertex_cover
  388. >>> print(list(independent_set))
  389. [2, 3, 4]
  390. See :mod:`bipartite documentation <networkx.algorithms.bipartite>`
  391. for further details on how bipartite graphs are handled in NetworkX.
  392. """
  393. # This is a Python implementation of the algorithm described at
  394. # <https://en.wikipedia.org/wiki/K%C3%B6nig%27s_theorem_%28graph_theory%29#Proof>.
  395. L, R = bipartite_sets(G, top_nodes)
  396. # Let U be the set of unmatched vertices in the left vertex set.
  397. unmatched_vertices = set(G) - set(matching)
  398. U = unmatched_vertices & L
  399. # Let Z be the set of vertices that are either in U or are connected to U
  400. # by alternating paths.
  401. Z = _connected_by_alternating_paths(G, matching, U)
  402. # At this point, every edge either has a right endpoint in Z or a left
  403. # endpoint not in Z. This gives us the vertex cover.
  404. return (L - Z) | (R & Z)
  405. #: Returns the maximum cardinality matching in the given bipartite graph.
  406. #:
  407. #: This function is simply an alias for :func:`hopcroft_karp_matching`.
  408. maximum_matching = hopcroft_karp_matching
  409. def minimum_weight_full_matching(G, top_nodes=None, weight="weight"):
  410. r"""Returns a minimum weight full matching of the bipartite graph `G`.
  411. Let :math:`G = ((U, V), E)` be a weighted bipartite graph with real weights
  412. :math:`w : E \to \mathbb{R}`. This function then produces a matching
  413. :math:`M \subseteq E` with cardinality
  414. .. math::
  415. \lvert M \rvert = \min(\lvert U \rvert, \lvert V \rvert),
  416. which minimizes the sum of the weights of the edges included in the
  417. matching, :math:`\sum_{e \in M} w(e)`, or raises an error if no such
  418. matching exists.
  419. When :math:`\lvert U \rvert = \lvert V \rvert`, this is commonly
  420. referred to as a perfect matching; here, since we allow
  421. :math:`\lvert U \rvert` and :math:`\lvert V \rvert` to differ, we
  422. follow Karp [1]_ and refer to the matching as *full*.
  423. Parameters
  424. ----------
  425. G : NetworkX graph
  426. Undirected bipartite graph
  427. top_nodes : container
  428. Container with all nodes in one bipartite node set. If not supplied
  429. it will be computed.
  430. weight : string, optional (default='weight')
  431. The edge data key used to provide each value in the matrix.
  432. Returns
  433. -------
  434. matches : dictionary
  435. The matching is returned as a dictionary, `matches`, such that
  436. ``matches[v] == w`` if node `v` is matched to node `w`. Unmatched
  437. nodes do not occur as a key in `matches`.
  438. Raises
  439. ------
  440. ValueError
  441. Raised if no full matching exists.
  442. ImportError
  443. Raised if SciPy is not available.
  444. Notes
  445. -----
  446. The problem of determining a minimum weight full matching is also known as
  447. the rectangular linear assignment problem. This implementation defers the
  448. calculation of the assignment to SciPy.
  449. References
  450. ----------
  451. .. [1] Richard Manning Karp:
  452. An algorithm to Solve the m x n Assignment Problem in Expected Time
  453. O(mn log n).
  454. Networks, 10(2):143–152, 1980.
  455. """
  456. import numpy as np
  457. import scipy as sp
  458. import scipy.optimize # call as sp.optimize
  459. left, right = nx.bipartite.sets(G, top_nodes)
  460. U = list(left)
  461. V = list(right)
  462. # We explicitly create the biadjancency matrix having infinities
  463. # where edges are missing (as opposed to zeros, which is what one would
  464. # get by using toarray on the sparse matrix).
  465. weights_sparse = biadjacency_matrix(
  466. G, row_order=U, column_order=V, weight=weight, format="coo"
  467. )
  468. weights = np.full(weights_sparse.shape, np.inf)
  469. weights[weights_sparse.row, weights_sparse.col] = weights_sparse.data
  470. left_matches = sp.optimize.linear_sum_assignment(weights)
  471. d = {U[u]: V[v] for u, v in zip(*left_matches)}
  472. # d will contain the matching from edges in left to right; we need to
  473. # add the ones from right to left as well.
  474. d.update({v: u for u, v in d.items()})
  475. return d