Construct adjacency matrix in MATLAB

Question

Consider a set of points arranged on a grid of size N-by-M. I am trying to build the adjacency matrix such that neighboring points are connected.

For example, in a 3x3 grid with a graph:

1-2-3
| | |
4-5-6
| | |
7-8-9

We should have the corresponding adjacency matrix:

+---+------------------------------------------------------+
|   |   1     2     3     4     5     6     7     8     9  |
+---+------------------------------------------------------+
| 1 |   0     1     0     1     0     0     0     0     0  |
| 2 |   1     0     1     0     1     0     0     0     0  |
| 3 |   0     1     0     0     0     1     0     0     0  |
| 4 |   1     0     0     0     1     0     1     0     0  |
| 5 |   0     1     0     1     0     1     0     1     0  |
| 6 |   0     0     1     0     1     0     0     0     1  |
| 7 |   0     0     0     1     0     0     0     1     0  |
| 8 |   0     0     0     0     1     0     1     0     1  |
| 9 |   0     0     0     0     0     1     0     1     0  |
+---+------------------------------------------------------+

As a bonus, the solution should work for both 4- and 8-connected neighboring points, that is:

   o             o  o  o
o  X  o   vs.    o  X  o
   o             o  o  o

This the code that I have so far:

N = 3; M = 3;
adj = zeros(N*M);

for i=1:N
    for j=1:M
        k = sub2ind([N M],i,j);
        if i>1
            ii=i-1; jj=j;
            adj(k,sub2ind([N M],ii,jj)) = 1; 
        end
        if i<N
            ii=i+1; jj=j;
            adj(k,sub2ind([N M],ii,jj)) = 1; 
        end
        if j>1
            ii=i; jj=j-1;
            adj(k,sub2ind([N M],ii,jj)) = 1; 
        end
        if j<M
            ii=i; jj=j+1;
            adj(k,sub2ind([N M],ii,jj)) = 1; 
        end
    end
end

How can this improved to avoid all the looping?

Answer 1

If you notice, there is a distinct pattern to the adjacency matrices you are creating. Specifically, they are symmetric and banded. You can take advantage of this fact to easily create your matrices using the diag function (or the spdiags function if you want to make a sparse matrix). Here is how you can create the adjacency matrix for each case, using your sample matrix above as an example:

4-connected neighbors:

mat = [1 2 3; 4 5 6; 7 8 9];                 % Sample matrix
[r, c] = size(mat);                          % Get the matrix size
diagVec1 = repmat([ones(c-1, 1); 0], r, 1);  % Make the first diagonal vector
                                             %   (for horizontal connections)
diagVec1 = diagVec1(1:end-1);                % Remove the last value
diagVec2 = ones(c*(r-1), 1);                 % Make the second diagonal vector
                                             %   (for vertical connections)
adj = diag(diagVec1, 1)+diag(diagVec2, c);   % Add the diagonals to a zero matrix
adj = adj+adj.';                             % Add the matrix to a transposed copy of
                                             %   itself to make it symmetric

And you'll get the following matrix:

adj =

     0  1  0  1  0  0  0  0  0
     1  0  1  0  1  0  0  0  0
     0  1  0  0  0  1  0  0  0
     1  0  0  0  1  0  1  0  0
     0  1  0  1  0  1  0  1  0
     0  0  1  0  1  0  0  0  1
     0  0  0  1  0  0  0  1  0
     0  0  0  0  1  0  1  0  1
     0  0  0  0  0  1  0  1  0

8-connected neighbors:

mat = [1 2 3; 4 5 6; 7 8 9];                 % Sample matrix
[r, c] = size(mat);                          % Get the matrix size
diagVec1 = repmat([ones(c-1, 1); 0], r, 1);  % Make the first diagonal vector
                                             %   (for horizontal connections)
diagVec1 = diagVec1(1:end-1);                % Remove the last value
diagVec2 = [0; diagVec1(1:(c*(r-1)))];       % Make the second diagonal vector
                                             %   (for anti-diagonal connections)
diagVec3 = ones(c*(r-1), 1);                 % Make the third diagonal vector
                                             %   (for vertical connections)
diagVec4 = diagVec2(2:end-1);                % Make the fourth diagonal vector
                                             %   (for diagonal connections)
adj = diag(diagVec1, 1)+...                  % Add the diagonals to a zero matrix
      diag(diagVec2, c-1)+...
      diag(diagVec3, c)+...
      diag(diagVec4, c+1);
adj = adj+adj.';                             % Add the matrix to a transposed copy of
                                             %   itself to make it symmetric

And you'll get the following matrix:

adj =

     0  1  0  1  1  0  0  0  0
     1  0  1  1  1  1  0  0  0
     0  1  0  0  1  1  0  0  0
     1  1  0  0  1  0  1  1  0
     1  1  1  1  0  1  1  1  1
     0  1  1  0  1  0  0  1  1
     0  0  0  1  1  0  0  1  0
     0  0  0  1  1  1  1  0  1
     0  0  0  0  1  1  0  1  0

Answer 2

Just for fun, here's a solution to construct the adjacency matrix by computing the distance between all pairs of points on the grid (not the most efficient way obviously)

N = 3; M = 3;                  %# grid size
CONNECTED = 8;                 %# 4-/8- connected points

%# which distance function
if CONNECTED == 4,     distFunc = 'cityblock';
elseif CONNECTED == 8, distFunc = 'chebychev'; end

%# compute adjacency matrix
[X Y] = meshgrid(1:N,1:M);
X = X(:); Y = Y(:);
adj = squareform( pdist([X Y], distFunc) == 1 );

And here's some code to visualize the adjacency matrix and the graph of connected points:

%# plot adjacency matrix
subplot(121), spy(adj)

%# plot connected points on grid
[xx yy] = gplot(adj, [X Y]);
subplot(122), plot(xx, yy, 'ks-', 'MarkerFaceColor','r')
axis([0 N+1 0 M+1])
%# add labels
[X Y] = meshgrid(1:N,1:M);
X = reshape(X',[],1) + 0.1; Y = reshape(Y',[],1) + 0.1;
text(X, Y(end:-1:1), cellstr(num2str((1:N*M)')) )

Answer 3

I just found this question when searching for the same problem. However, none of the provided solutions worked for me because of the problem size which required the use of sparse matrix types. Here is my solution which works on large scale instances:

function W = getAdjacencyMatrix(I)

[m, n] = size(I);

I_size = m*n;

% 1-off diagonal elements
V = repmat([ones(m-1,1); 0],n, 1);
V = V(1:end-1); % remove last zero

% n-off diagonal elements
U = ones(m*(n-1), 1);

% get the upper triangular part of the matrix
W = sparse(1:(I_size-1),    2:I_size, V, I_size, I_size)...
  + sparse(1:(I_size-m),(m+1):I_size, U, I_size, I_size);

% finally make W symmetric
W = W + W';

Answer 4

Just came across this question. I have a nice working m-function (link: sparse_adj_matrix.m) that is quite general.

It can handle 4-connect grid (radius 1 according to L1 norm), 8-connect grid (radius 1 according to L_infty norm).
It can also support 3D (and arbitrarily higher domensional grids).
The function can also connect nodes further than radius = 1.

Here's the signiture of the function:

% Construct sparse adjacency matrix (provides ii and jj indices into the
% matrix)
%
% Usage:
%   [ii jj] = sparse_adj_matrix(sz, r, p)
%
% inputs:
%   sz - grid size (determine the number of variables n=prod(sz), and the
%        geometry/dimensionality)
%   r  - the radius around each point for which edges are formed
%   p  - in what p-norm to measure the r-ball, can be 1,2 or 'inf'
%
% outputs
%   ii, jj - linear indices into adjacency matrix (for each pair (m,n)
%   there is also the pair (n,m))
%
% How to construct the adjacency matrix?
% >> A = sparse(ii, jj, ones(1,numel(ii)), prod(sz), prod(sz));
%
%
% Example:
% >> [ii jj] = sparse_adj_matrix([10 20], 1, inf);
% construct indices for 200x200 adjacency matrix for 8-connect graph over a
% grid of 10x20 nodes.
% To visualize the graph:
% >> [r c]=ndgrid(1:10,1:20);
% >> A = sparse(ii, jj, 1, 200, 200);;
% >> gplot(A, [r(:) c(:)]);

Answer 5

Your current code doesn't seem so bad. One way or another you need to iterate over all neighbor pairs. If you really need to optimize the code, I would suggest:

• loop over node indices i, where 1 <= i <= (N*M)
• don't use sub2ind() for efficiency, the neighbors of node i are simpy [i-M, i+1, i+M, i-1] in clockwise order

Notice that to get all neighbor pairs of nodes:

• you only have to compute the "right" neighbors (i.e. horizontal edges) for nodes i % M != 0 (since Matlab isn't 0-based but 1-based)
• you only have to compute "above" neighbors (i.e. vertical edges) for nodes i > M
• there is a similar rule for diagonal edges

This would leed to a single loop (but same number of N*M iterations), doesn't call sub2ind(), and has only two if statements in the loop.

Answer 6

For each node in the graph add a connection to the right and one downwards. Check that you don't overreach your grid. Consider the following function that builds the adjacency matrix.

function  adj = AdjMatrixLattice4( N, M )
    % Size of adjacency matrix
    MN = M*N;
    adj = zeros(MN,MN);

    % number nodes as such
    %  [1]---[2]-- .. --[M]
    %   |     |          |
    % [M+1]-[M+2]- .. -[2*M]
    %   :     :          :
    %   []    []   ..  [M*N]     

    for i=1:N
        for j=1:N
            A = M*(i-1)+j;          %Node # for (i,j) node
            if(j<N)                
                B = M*(i-1)+j+1;    %Node # for node to the right
                adj(A,B) = 1;
                adj(B,A) = 1;
            end
            if(i<M)
                B = M*i+j;          %Node # for node below
                adj(A,B) = 1;       
                adj(B,A) = 1;
            end            
        end
    end    
end

Example as above AdjMatrixLattice4(3,3)=

 0     1     0     1     0     0     0     0     0
 1     0     1     0     1     0     0     0     0
 0     1     0     0     0     1     0     0     0
 1     0     0     0     1     0     1     0     0
 0     1     0     1     0     1     0     1     0
 0     0     1     0     1     0     0     0     1
 0     0     0     1     0     0     0     1     0
 0     0     0     0     1     0     1     0     1
 0     0     0     0     0     1     0     1     0