在 MATLAB 中构造邻接矩阵

Question

考虑排列在大小为 N×M 的网格上的一组点。 我正在尝试构建邻接矩阵 相邻点相连。

例如，在带有图的 3x3 网格中：

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

我们应该有相应的邻接矩阵：

+---+------------------------------------------------------+
|   |   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  |
+---+------------------------------------------------------+

作为奖励，该解决方案应该适用于 4 个和 8 个连接的相邻点，即：

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

这是我到目前为止的代码：

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

如何改进以避免所有循环？

Answer 1

如果您注意到，您正在创建的邻接矩阵有一个独特的模式。具体来说，它们是对称的并且带状。您可以利用这一事实，使用 diag 轻松创建矩阵 函数（或 spdiags函数（如果你想创建一个稀疏矩阵）。以下是如何为每种情况创建邻接矩阵，以上面的示例矩阵为例：

4 个连接的邻居：

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

您将得到以下矩阵：

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 个连接的邻居：

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

您将得到以下矩阵：

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

只是为了好玩，这里有一个通过计算网格上所有点对之间的距离来构造邻接矩阵的解决方案（显然不是最有效的方法）

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 );

这里有一些代码来可视化邻接矩阵和连接点的图：

%# 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

我刚刚在搜索同样的问题时发现了这个问题。然而，由于问题大小需要使用稀疏矩阵类型，所提供的解决方案都不适合我。这是我的解决方案，适用于大规模实例：

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

刚刚遇到这个问题。我有一个很好的工作 m 函数（链接：sparse_adj_matrix.m）这是非常普遍的。

它可以处理4连接网格（根据L1规范半径1），8连接网格（根据L_infty规范半径1）。
它还可以支持 3D（以及任意更高的维度网格）。
该函数还可以连接比 radius = 1 更远的节点。

这是该函数的签名：

% 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

您当前的代码看起来并没有那么糟糕。您需要以某种方式迭代所有邻居对。如果您确实需要优化代码，我建议：

• 循环遍历节点索引 i，其中 1 <= i <= (N*M)
• 不要使用 sub2ind() 来提高效率，节点 i 的邻居按顺时针顺序简单地 [iM, i+1, i+M, i-1]

请注意，要获取节点的所有邻居对：

• 您只需计算“右” " 节点 i % M != 0 的邻居（即水平边缘）（因为 Matlab 不是基于 0 而是基于 1），
• 您只需计算“上方”邻居（即垂直边缘）对于节点 i > M
• 对于对角线边缘有类似的规则

这将导致单个循环（但 N*M 迭代次数相同），不会调用 sub2ind()，并且循环中只有两个 if 语句。

Answer 6

好吧，像往常一样，我迟到了，但是朋友之间的十年又算什么呢？
再加上这个答案基于（吹喇叭）RECURSION，并将上述问题推广到大小为 N1 x N2 x ... x Nd 的任意 d 维晶格。
该解决方案本身就是 https://stackoverflow.com/a/77100541/4474463
我倾向于以教学方式写作，所以不要反对我，请原谅这个符号，因为这是少数不支持 LaTeX 的 SO 之一。
我们首先研究一个简单的例子并用图片来思考。
然后我们将这些图片转化为数学，最后将数学转化为算法。
准备好？

想象一个 4x3x2 的格子。我们如何才能“从头开始”构建它？
一种方法是从一维晶格开始：具有四个顶点和 3 条边的图：

酷，这很容易。现在让我们将其扩展到 4x3 的二维晶格。
为此，我们制作 1d 晶格的 3-1 个副本，平移它们，然后将所有副本粘合在一起形成 2d 晶格：

我打赌你能明白这是怎么回事。最后，我们制作 2-1 个 2d 晶格的副本，并将它们粘合在一起：

哒哒！

应该清楚的是，您可以类似地制作任意尺寸的 d 维晶格。
现在，我们如何将这个图形过程转化为数学？
我们需要工具来复制邻接矩阵并将其粘合在一起。
对我们来说幸运的是，我们所需要的只是 Kronecker 产品，它具有实现这两点的强大功能。
而且它已经在每个值得一试的线性代数库中实现了。

在向您展示如何使用克罗内克积之前，我们需要首先了解长度为 N 的一维晶格的邻接矩阵。
这应该不会太难：顶点 1 仅连接到顶点 2，顶点 N 仅连接到 N-1，顶点 k 之间的每个顶点都连接到 k+1 和 k-1。
所以邻接矩阵为：

另外，对角线上有 N 个 1 的单位矩阵表示为 1(N)，并且给定 M(k-1) = N1 * N2 * ... * N(k-1)。

假设您知道 N1 x N2 x ... x N(k-1) 格的邻接矩阵，将其称为 A(k-1)。
然后，为了制作 N1 x N2 x ... x N(k-1) x Nk 格的邻接矩阵，我们首先采用 A(k-1) 的 Nk 个副本。
然后我们在每个副本中取相同的顶点，并将它们与边顺序粘合在一起以形成一维晶格。
生成的邻接矩阵具有以下形式

如果您在解析右侧时遇到困难，第一项对应于 (k-1)d 晶格的 Nk 个副本；这就是中间矩阵中的所有对角项。
请注意所有 1(M(k-1)) 的位置与我们刚刚看到的一维晶格的位置相同：这些是将复制的顶点粘合在一起的边。
就是这样，我们有一种构造格子的递归方法。
当然，基本情况是一维晶格。

这是 MATLAB 实现，您可以将其翻译成您最喜欢的编程语言。它以正整数、dim（每个维度中的格点数）组成的向量作为输入，并吐出邻接矩阵。

function L = LATTICE(dims)
  if length(dims) == 1
     L = sparse(eye(dims));
     L = circshift(L,1,1) + circshift(L,1,2);
     if dims > 2
        L(1,end) = 0; L(end,1) = 0;
     end
  else
     m = prod(dims(1:(end-1));
     L = kron(eye(dims(end)),LATTICE(dims(1:(end-1)) + ...
             kron(LATTICE(dims(end)),eye(m));
  end
end

注意稀疏矩阵的使用。你应该使用任何类似的结构，只是不要制作一个巨大的格子并尝试查看完整的矩阵，好吗？

使用上面的内容制作 12 x 8 x 6 格子的一个例子是

L = LATTICE([12,8,6])

，绘制它

plot(graph(L),"Layout", "force3")

给我们这个漂亮的野兽：

（当然，经过一些表面上的调整。）

希望您发现这个回复既有用又读起来令人愉快。
干杯!

Answer 7

对于图中的每个节点，添加一个向右的连接和一个向下的连接。检查并确保您没有超出网格范围。考虑以下构建邻接矩阵的函数。

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

如上示例 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