function [l v] = grasshopper(W, r, lambda, k)
%
% Reranking items by random walk on graph with absorbing states.
%
% INPUT
%
% W: n*n weight matrix with non-negative entries.
% W(i,j) = edge weight for the edge i->j
% The graph can be directed or undirected. Self-edges are allowed.
% r: user-supplied initial ranking of the n items.
% r is a probability vector of size n.
% highly ranked items have larger r values.
% lambda:
% trade-off between relying completely on W and ignoring r (lambda=1)
% and completely relying on r and ignoring W (lambda=0).
% k: return top k items as ranked by Grasshopper.
%
% OUTPUT
%
% l: the top k indices after reranking.
% l(1) is the best item's index (into W), l(2) the second best, and so on
% v: v(1) is the first item's stationary probability;
% v(2)..v(k) are the next k-1 item's average number of visits during
% absorbing random walks, in the respective iterations when they were
% selected.
n = size(W,1);
% sanity check first
if (min(min(W))<0),
error('Found a negative entry in W! Stop.');
elseif (abs(sum(r)-1)>1e-5),
error('The vector r does not sum to 1! Stop.');
elseif (lambda<0 | lambda>1)
error('lambda is not in [0,1]! Stop.');
elseif (k<0 | k>n)
error('k is not in [0,n]! Stop.');
end
% creating the graph-based transition matrix P
P = W ./ repmat(sum(W,2), 1, n); % row normalized
% incorporating user-supplied initial ranking by adding teleporting
hatP = lambda*P + (1-lambda)*repmat(r', n, 1);
% finding the highest ranked item from the stationary distribution
% of hatP: q=hatP'*q. The item with the largest stationary probability
% is our higest ranked item.
q = stationary(hatP);
% the most probably state is our first absorbing node. Put it in l.
[v, l]=max(q);
% reranking k-1 more items by picking out the most-visited node one by one.
% once picked out, the item turns into an absorbing node.
while (length(l)