|
| 1 | +from collections import OrderedDict |
| 2 | +from warnings import warn |
| 3 | + |
| 4 | +import grblas as gb |
| 5 | +import networkx as nx |
| 6 | +from grblas import Vector, binary, unary |
| 7 | +from grblas.semiring import plus_first, plus_times |
| 8 | + |
| 9 | + |
| 10 | +def pagerank_core( |
| 11 | + A, |
| 12 | + alpha=0.85, |
| 13 | + personalization=None, |
| 14 | + max_iter=100, |
| 15 | + tol=1e-06, |
| 16 | + nstart=None, |
| 17 | + dangling=None, |
| 18 | + row_degrees=None, |
| 19 | + name="pagerank", |
| 20 | +): |
| 21 | + N = A.nrows |
| 22 | + if A.nvals == 0: |
| 23 | + return Vector.new(float, N, name=name) |
| 24 | + |
| 25 | + # Initial vector |
| 26 | + x = Vector.new(float, N, name="x") |
| 27 | + if nstart is None: |
| 28 | + x[:] = 1.0 / N |
| 29 | + else: |
| 30 | + denom = nstart.reduce(allow_empty=False).value |
| 31 | + if denom == 0: |
| 32 | + raise ZeroDivisionError() |
| 33 | + x << nstart / denom |
| 34 | + |
| 35 | + # Personalization vector or scalar |
| 36 | + if personalization is None: |
| 37 | + p = 1.0 / N |
| 38 | + else: |
| 39 | + denom = personalization.reduce(allow_empty=False).value |
| 40 | + if denom == 0: |
| 41 | + raise ZeroDivisionError() |
| 42 | + p = (personalization / denom).new(name="p") |
| 43 | + |
| 44 | + # Inverse of row_degrees |
| 45 | + # Fold alpha constant into S |
| 46 | + if row_degrees is None: |
| 47 | + S = A.reduce_rowwise().new(float, name="S") |
| 48 | + S << alpha / S |
| 49 | + else: |
| 50 | + S = (alpha / row_degrees).new(name="S") |
| 51 | + |
| 52 | + if A.ss.is_iso: |
| 53 | + # Fold iso-value of A into S |
| 54 | + # This lets us use the plus_first semiring, which is faster |
| 55 | + iso_value = A.ss.iso_value |
| 56 | + if iso_value != 1: |
| 57 | + S *= iso_value |
| 58 | + semiring = plus_first[float] |
| 59 | + else: |
| 60 | + semiring = plus_times[float] |
| 61 | + |
| 62 | + is_dangling = S.nvals < N |
| 63 | + if is_dangling: |
| 64 | + dangling_mask = Vector.new(float, N, name="dangling_mask") |
| 65 | + dangling_mask(mask=~S.S) << 1.0 |
| 66 | + # Fold alpha constant into dangling_weights (or dangling_mask) |
| 67 | + if dangling is not None: |
| 68 | + dangling_weights = (alpha / dangling.reduce(allow_empty=False).value * dangling).new( |
| 69 | + name="dangling_weights" |
| 70 | + ) |
| 71 | + elif personalization is None: |
| 72 | + # Fast case (and common case); is iso-valued |
| 73 | + dangling_mask(mask=dangling_mask.S) << alpha * p |
| 74 | + else: |
| 75 | + dangling_weights = (alpha * p).new(name="dangling_weights") |
| 76 | + |
| 77 | + # Fold constant into p |
| 78 | + p *= 1 - alpha |
| 79 | + |
| 80 | + # Power iteration: make up to max_iter iterations |
| 81 | + xprev = Vector.new(float, N, name="x_prev") |
| 82 | + w = Vector.new(float, N, name="w") |
| 83 | + for _ in range(max_iter): |
| 84 | + xprev, x = x, xprev |
| 85 | + |
| 86 | + # x << alpha * ((xprev * S) @ A + "dangling_weights") + (1 - alpha) * p |
| 87 | + x << p |
| 88 | + if is_dangling: |
| 89 | + if dangling is None and personalization is None: |
| 90 | + # Fast case: add a scalar; x is still iso-valued (b/c p is also scalar) |
| 91 | + x += xprev @ dangling_mask |
| 92 | + else: |
| 93 | + # Add a vector |
| 94 | + x += plus_first(xprev @ dangling_mask) * dangling_weights |
| 95 | + w << xprev * S |
| 96 | + x += semiring(w @ A) # plus_first if A.ss.is_iso else plus_times |
| 97 | + |
| 98 | + # Check convergence, l1 norm: err = sum(abs(xprev - x)) |
| 99 | + xprev << binary.minus(xprev | x, require_monoid=False) |
| 100 | + xprev << unary.abs(xprev) |
| 101 | + err = xprev.reduce().value |
| 102 | + if err < N * tol: |
| 103 | + x.name = name |
| 104 | + return x |
| 105 | + raise nx.PowerIterationFailedConvergence(max_iter) |
| 106 | + |
| 107 | + |
| 108 | +def pagerank( |
| 109 | + G, |
| 110 | + alpha=0.85, |
| 111 | + personalization=None, |
| 112 | + max_iter=100, |
| 113 | + tol=1e-06, |
| 114 | + nstart=None, |
| 115 | + weight="weight", |
| 116 | + dangling=None, |
| 117 | +): |
| 118 | + warn("", DeprecationWarning, stacklevel=2) |
| 119 | + N = len(G) |
| 120 | + if N == 0: |
| 121 | + return {} |
| 122 | + node_ids = OrderedDict((k, i) for i, k in enumerate(G)) |
| 123 | + A = gb.io.from_networkx(G, nodelist=node_ids, weight=weight, dtype=float) |
| 124 | + |
| 125 | + x = p = dangling_weights = None |
| 126 | + # Initial vector (we'll normalize later) |
| 127 | + if nstart is not None: |
| 128 | + indices, values = zip(*((node_ids[key], val) for key, val in nstart.items())) |
| 129 | + x = Vector.from_values(indices, values, size=N, dtype=float, name="nstart") |
| 130 | + # Personalization vector (we'll normalize later) |
| 131 | + if personalization is not None: |
| 132 | + indices, values = zip(*((node_ids[key], val) for key, val in personalization.items())) |
| 133 | + p = Vector.from_values(indices, values, size=N, dtype=float, name="personalization") |
| 134 | + # Dangling nodes (we'll normalize later) |
| 135 | + row_degrees = A.reduce_rowwise().new(name="row_degrees") |
| 136 | + if dangling is not None: |
| 137 | + if row_degrees.nvals < N: # is_dangling |
| 138 | + indices, values = zip(*((node_ids[key], val) for key, val in dangling.items())) |
| 139 | + dangling_weights = Vector.from_values( |
| 140 | + indices, values, size=N, dtype=float, name="dangling" |
| 141 | + ) |
| 142 | + result = pagerank_core( |
| 143 | + A, |
| 144 | + alpha=alpha, |
| 145 | + personalization=p, |
| 146 | + max_iter=max_iter, |
| 147 | + tol=tol, |
| 148 | + nstart=x, |
| 149 | + dangling=dangling_weights, |
| 150 | + row_degrees=row_degrees, |
| 151 | + ) |
| 152 | + if result.nvals != N: |
| 153 | + # Not likely, but fill with 0 just in case |
| 154 | + result(mask=~result.S) << 0 |
| 155 | + return dict(zip(node_ids, result.to_values()[1])) |
0 commit comments