Hits.m

function [ auth, hub ] = HITS(a, epsilon)
%HITS Computes HITS, returned as a numeric array.
%   a is an adjacent sparse matrix representing a graph
%   epsilone defines the threshold to stop

n = size(a, 1);

% Initialize authority scores and hub scores
auth = ones(n, 1);
hub = ones(n, 1);

% Initialize error threshold
delta_auth = Inf;
delta_hub = Inf;
while delta_auth > epsilon && delta_hub > epsilon
    % Calc authority scores
    auth_prev = auth;
    auth = a.' * (a * auth);
    auth = auth / norm(auth,1); % Normalize
    % Compute error for authority
    delta_auth = norm(auth - auth_prev, 1);
    
    % Calc hub score
    hub_prev = hub;
    hub = a * (a.' * hub);
    hub = hub / norm(hub,1); % Normalize
    
    % Compute error for hub
    delta_hub = norm(hub - hub_prev, 1);
end
end