First version of Anymatrix

co-authored-by: Nicholas. J. Higham <[email protected]>
co-authored-by: Mantas Mikaitis <[email protected]>

Author: Mantas Mikaitis

anymatrix.m

@@ -0,0 +1,6 @@
+% Template for an alias of anymatrix.
+% Rename the file and the function name as appropriate.
+% For example am.m and am().
+function varargout = anymatrix_alias(varargin)
+[varargout{1:nargout}] = anymatrix(varargin{1:nargin});
+end

anymatrix_alias.m

@@ -0,0 +1,6 @@
+% A function that links anymatrix to the matrix-generating functions
+% inside the /private folders of each group.
+function varargout = anymatrix_contest(matrix_name, varargin)
+handle = str2func(matrix_name);
+[varargout{1:nargout}] = handle(varargin{1:nargin-1});
+end

contest/anymatrix_contest.m

@@ -0,0 +1,28 @@
+% CONTEST.
+% Version 1.2                20-May-2008
+% Alan Taylor and Des Higham
+%
+% baitsample    -   Bait and prey subsampling
+% curvature     -   Curvatures (clustering coefficients)
+% erdrey        -   Erdos-Renyi model graph
+% geo           -   Geometric random graph
+% gilbert       -   Gilbert model graph
+% kleinberg     -   Kleinberg model graph
+% lap           -   Laplacian matrix (normalized or unnormalized)
+% lockandkey    -   Lock and key model graph
+% mht           -   Mean hitting times
+% pagerank      -   PageRank matrix
+% pathlength    -   Minimum path lengths
+% pref          -   Scale free random graph
+% renga         -   Range dependent random graph
+% rewire        -   Redirect edges
+% short         -   Add shortcuts
+% smallw        -   Small world random graph
+% sticky        -   Stickiness model random graph
+% unisample     -   Uniform subsampling
+For documentation see 
+- https://www.maths.ed.ac.uk/~dhigham/CONTEST_package.html
+- Alan Taylor and Desmond J. Higham. CONTEST: A Controllable Test Matrix
+  Toolbox for MATLAB. ACM Trans. Math. Software, 35(4):26:1--26:17, 2009.
+  https://doi.org/10.1145/1462173.1462175
+Included with permission.

contest/private/Contents.m

@@ -0,0 +1,22 @@
+function P = am_properties
+
+P = {
+    'baitsample', {'random', 'symmetric', 'sparse'}
+    'curvature', {}
+    'erdrey', {'random', 'symmetric', 'sparse'}
+    'geo', {'random', 'symmetric', 'sparse'}
+    'gilbert', {'random', 'symmetric', 'sparse'}
+    'kleinberg', {'random', 'symmetric', 'sparse'}
+    'lap', {'random', 'symmetric', 'sparse'}
+    'lockandkey', {'random', 'symmetric', 'sparse'}
+    'mht', {'random', 'sparse'}
+    'pagerank', {'random', 'sparse'}
+    'pathlength', {'random', 'symmetric', 'sparse'}
+    'pref', {'random', 'symmetric', 'sparse'}
+    'renga', {'random', 'symmetric', 'sparse'}
+    'rewire', {'random', 'symmetric', 'sparse'}
+    'short', {'random', 'symmetric', 'sparse'}
+    'smallw', {'random', 'symmetric', 'sparse'}
+    'sticky', {'random', 'symmetric', 'sparse'}
+    'unisample', {'random', 'symmetric', 'sparse'}
+};

contest/private/am_properties.m

@@ -0,0 +1,67 @@
+function B = baitsample(A,bait,prey)
+
+%BAITSAMPLE     bait and prey subsampling
+%
+%   Input   A: n by n adjacency matrix
+%           bait: proportion of nodes to sample. Defaults to 0.5.
+%           prey: proportion of edges to retain from each "bait" node.
+%                 Defaults to 0.5.
+%
+%   Output  B: adjacency matrix with the attribute sparse.
+%              Dimension of B cannot be predicted.
+%
+%   Description:    The adjacency matrix A is considered, and a 
+%                   proportion, bait, of its rows/columns are retained.
+%                   Of these, a proportion, prey, of the outgoing edges
+%                   are retained. All other entries are set to zero.
+%                   Disconnected nodes are then removed.
+%
+%   Reference:  J. Han, D. Dupuy, N. Bertin, M. Cusick, M. Vidal,
+%               Effect of sampling on topology predictions of
+%               protein-protein interaction networks,
+%               Nature Biotechnology 23 (2005), pp. 839-844.
+%
+%   Example: B = baitsample(A,0.1,0.4);
+
+if nargin <= 2
+    prey = 0.5;
+    if nargin == 1
+        bait = 0.5;
+    end
+end
+
+n = length(A);
+B = sparse(n,n);
+
+rp = randperm(n);
+ibait = rp(1:ceil(bait*n));     % bait rows/columns to be retained
+
+for i = 1:length(ibait)
+    
+    B(ibait(i),:) = A(ibait(i),:);
+    B(:,ibait(i)) = A(:,ibait(i));
+    iprey = find(B(ibait(i),:));
+    psum = ceil(prey*length(iprey)); % number of prey nodes for ith bait
+    
+    if length(iprey) ~= 0
+        dist = (1/length(iprey)) : (1/length(iprey)) : 1;
+    
+        while length(find(B(ibait(i),:))) > psum
+        
+            r = rand;
+            pos = min(find(r<=dist));
+            B(ibait(i),iprey(pos))=0;
+            B(iprey(pos),ibait(i))=0;
+            
+        end
+        
+    end
+    
+end
+
+for i = 1 : n               % trim isolated nodes
+	if ~any(B(:,(n-i+1)))
+		B(:,(n-i+1)) = [];
+		B((n-i+1),:) = [];
+	end
+end

contest/private/baitsample.m

@@ -0,0 +1,54 @@
+function C = curvature(A,ind)
+
+%CURVATURE  Compute the curvatures (clustering coefficients) for a given
+%           adjacency matrix.
+%
+%   Input   A:   n by n adjacency matrix.
+%           ind: node index (optional). Can also take string values
+%                'max' for maximum and 'ave' for average, ignoring those
+%                that are undefined.
+%
+%   Output  C: n by 1 vector of curvatures (scalar if ind is provided).
+%              Undefined values are returned as NaN.
+%
+%   Description:    Calculates the curvatures of the nodes in a graph.
+%                   Defined for each node by the frequency with which
+%                   its neighbours are themselves connected.
+%
+%   Examples: C = curvature(A)           vector of curvatures
+%             C = curvature(A,4)         curvature of 4th node
+%             C = curvature(A,'max')     largest curvature
+%             C = curvature(A,'ave')      average curvature
+
+n = length(A);
+
+v = sum(A);
+
+b = diag(A^3);
+d = v.*(v-1);
+
+ccoefs = zeros(1,n);
+
+for i = 1:n
+    if d(i) <= 1
+        ccoefs(i) = NaN;
+    else
+        ccoefs(i) = b(i)/d(i);
+    end
+end
+
+if nargin==1
+    C = ccoefs;
+else
+    if strcmp(ind,'max')
+        C = max(ccoefs);
+    else
+        
+        if strcmp(ind,'ave')
+            C = mean(ccoefs(find(~isnan(ccoefs))));
+        else
+        
+            C = ccoefs(ind);
+        end
+    end
+end

contest/private/curvature.m

@@ -0,0 +1,60 @@
+function A = erdrey(n,m)
+
+%ERDREY     Generate adjacency matrix for a G(n,m) type random graph.
+%
+%   Input   n: dimension of matrix (number of nodes in graph).
+%           m: 2*m is the number of 1's in matrix (number of edges in graph).
+%           Defaults to the smallest integer larger than n*log(n)/2.
+%
+%   Output  A: n by n symmetric matrix with the attribute sparse.
+%
+%
+%   Description:    An undirected graph is chosen uniformly at random from
+%                   the set of all symmetric graphs with n nodes and m
+%                   edges.
+%  
+%   Reference:  P. Erdos, A. Renyi,
+%               On Random Graphs,
+%               Publ. Math. Debrecen, 6 1959, pp. 290-297.
+%
+%   Example: A = erdrey(100,10);
+
+if nargin == 1
+    m = ceil(n*log(n)/2);
+end
+
+nonzeros = ceil(0.5*n*(n-1)*rand(m,1));
+v = zeros(n,1);
+for count = 1:n
+    v(count) = count*(count-1)/2;
+end
+
+I = zeros(m,1);
+J = zeros(m,1);
+S = ones(m,1);
+
+for count = 1:m
+    i = min(find(v >= nonzeros(count)));
+    j = nonzeros(count) - (i-1)*(i-2)/2;
+    I(count) = i;
+    J(count) = j;
+end
+
+A = sign(sparse([I;J],[J;I],[S;S],n,n));
+
+while nnz(A) ~= 2*m
+    
+    difference = m-nnz(A)/2;
+    Inew = zeros(difference,1);
+    Jnew = zeros(difference,1);
+    for count = 1:difference
+        index = ceil(0.5*n*(n-1)*rand);
+        Inew(count) = min(find(v>=index));
+        Jnew(count) = index - (Inew(count)-1)*(Inew(count)-2)/2;
+    end
+    I = cat(1,I,Inew);
+    J = cat(1,J,Jnew);
+    S = ones(length(I),1);
+    A = sign(sparse([I;J],[J;I],[S;S],n,n));
+    
+end

contest/private/erdrey.m

@@ -0,0 +1,71 @@
+function A = geo(n,r,m,per,pnorm)
+
+%GEO      Generate adjacency matrix for a geometric random graph.
+%
+%   Input   n: dimension of matrix (number of nodes in graph)
+%           r: radius used to defined entries (edges). Defaults to the
+%           square root of 1.44/n.
+%           m: dimension of coordinate system. Defaults to 2.
+%         per: periodicity of coordinate system. Periodic if per == 1,
+%              not periodic if per == 0. Defaults to 0.
+%       pnorm: norm to measure distance between nodes. Defaults to 2.
+%
+%   Output  A: n by n symmetric matrix with the attribute sparse
+%
+%
+%   Description: nodes are placed randomly in the unit m-cube.
+%                An edge is created if two nodes are within distance r.
+%   Reference:   M. Penrose, Geometric Random Graphs,
+%                Oxford Univeristy Press, 2003.
+%
+%   Example: A = geo(100,0.01,3,1,2);
+
+if nargin <= 4
+    pnorm = 2;
+    if nargin <= 3
+        per = 0;
+        if nargin <= 2
+            m = 2;
+            if nargin == 1
+                r = sqrt(1.44/n);
+            end
+        end
+    end
+end
+
+coords = rand(n,m);
+
+I = [];
+J = [];
+
+if per == 0
+    
+    for i = 2:n
+        for j = 1:(i-1)
+            diff = abs(coords(i,:) - coords(j,:));
+            if norm(diff,pnorm)<=r
+                J = cat(1,J,j);
+                I = cat(1,I,i);
+            end
+        end
+    end
+    
+end
+
+if per == 1
+    
+    for i = 2:n
+        for j = 1:(i-1)
+            diff = min( abs( coords(i,:) - coords(j,:) ), abs( 1 - abs( coords(i,:) - coords(j,:) ) ) );
+            if norm(diff,pnorm)<=r
+                J = cat(1,J,j);
+                I = cat(1,I,i);
+            end
+        end
+        
+    end
+    
+end
+
+S = ones(length(I),1);
+A = sparse([I;J],[J;I],[S;S],n,n);

contest/private/geo.m

@@ -0,0 +1,59 @@
+function A = gilbert(n,p);
+
+%GILBERT     Generate adjacency matrix for a G(n,p) type random graph.
+%
+%   Input    n: dimension of matrix (number of nodes in graph).
+%            p: probability that any two nodes are neighbours. Defaults to
+%               log(n)/n.
+%
+%   Output   A: n by n symmetric matrix with the attribute sparse.
+%
+%
+%   Description:    An undirected graph is created by considering pairs of
+%                   nodes and connecting them with independent probability p.
+%
+%   Reference:      E.N. Gilbert,
+%                   Random Graphs,
+%                   Ann. Math. Statist.,30, (1959) pp. 1141-1144.
+%
+%   This code is a direct translation of Algorithm 1 in
+%                   V. Batagelj, U. Brandes,
+%                   Efficient generation of large random networks,
+%                   Phys. Rev. E, 71 (2005).
+%   This algorithm uses a geometric method to skip over potential edges, 
+%   and has optimal complexity.
+%
+%   Example: A = gilbert(100,0.1);
+
+if nargin == 1
+    p = log(n)/n;
+end
+
+v = zeros(n,1);        % Think of lower triangle of n-by-n array ordered
+for k = 1:n,           % lexographically, row-wise. So v(k) is the biggest 
+    v(k) = k*(k-1)/2;  % index appearing in row k. 
+end
+
+I = zeros(ceil(0.5*p*n^2),1);  % We don't know exact length
+J = zeros(ceil(0.5*p*n^2),1);  % Expected length is 0.5*p*n(n-1)
+
+count = 0;
+w = 0;
+
+w = w + 1 + floor(log(1-rand)/log(1-p));  %Lexographical index of next nonzero
+while w < n*(n-1)/2;
+   i = min(find(v >= w));    % Recover i and j from
+   j = w - (i-1)*(i-2)/2;    % lexographical index w
+   I(count+1) = i;
+   J(count+1) = j;
+   count = count + 1;
+   w = w + 1 + floor(log(1-rand)/log(1-p)); %Lexographical index of next nonzero
+end
+
+Ifind = find(I>0);   % trim any left-over zeros 
+I = I(Ifind);        % from I and J
+Jfind = find(J>0);
+J = J(Jfind);
+   
+S = ones(length(I),1);
+A = sparse([I;J],[J;I],[S;S],n,n);