在MATLAB中构造3D图邻接矩阵

问题描述 投票:0回答:1

我对扩展这个问题/答案(https://stackoverflow.com/a/3283732/2371031)感兴趣,以将4个连接的案例扩展到第三维。

问题1:给定一个X x Y x Z维矩阵,如何构造4连接邻接矩阵?

问题2:给定邻接矩阵,如何计算边缘权重作为连接节点值的平均值?

问题3:如何解决非常大的矩阵的问题1和2而又不会耗尽内存?稀疏矩阵或边列表似乎是可能的解决方案,但如何实现?

% simple case
X = 2; Y = 2; Z = 2;
mat = reshape([1:(X*Y*Z)],X,Y,Z);
>> mat

mat(:,:,1) =

     1     3
     2     4


mat(:,:,2) =

     5     7
     6     8


% large case
X = 256; Y = 256; Z = 1000;
mat = reshape([1:(X*Y*Z)],X,Y,Z);

matlab matrix graph-theory graph-algorithm
1个回答
0
投票

答案1:下面的答案3中的循环效率更高一些。与其建立巨大的边缘和权重矩阵,不如直接在循环中迭代地建立图形。这可以很快完成64x64x10(秒),但是对于大型示例(256x256x1000),运行时间要长得多。制定一个更快的版本会很好。

X=256;Y=256;Z=1000;
% input matrix of 0-1 values 
mat = rand(X,Y,Z);
% matrix of indices for input matrix coordinates
imat = reshape([1:(X*Y*Z)],X,Y,Z);

% initialize empty graph
g = digraph();

xdim = size(imat,1);
ydim = size(imat,2);
zdim = size(imat,3);

% create edge list with weights
disp('calculating start/end nodes for edges and associated weights...')


for y = 1:ydim  % loop through all "horizontal" (column) connections in each z plane
    % example:
    % 1 - 3
    % 2 - 4
    for z = 1:zdim
        if y ~= 1
            nodes1 = imat(:,y,z);   % grab all rows (nodes) of column y 
            nodes2 = imat(:,y-1,z); % grab all rows (nodes) of column y-1
            wts = ((1-mat(:,y,z)) + (1-mat(:,y-1,z))) ./ 2;  % average of 1-value for each node
            g = addedge(g,nodes1,nodes2,wts); % add to graph
        end
        if y ~= ydim
            nodes1 = imat(:,y,z);
            nodes2 = imat(:,y+1,z);
            wts = ((1-mat(:,y,z)) + (1-mat(:,y+1,z))) ./ 2;  % average of 1-value for each node
            g = addedge(g,nodes1,nodes2,wts); % add to graph
        end
    end
end

for x = 1:xdim  % loop through all "vertical" (row) connections within each z plane
    % example:
    % 1  3
    % |  |
    % 2  4

    for z = 1:zdim
        if x ~= 1
            nodes1 = imat(x,:,z);   % grab all rows (nodes) of column y 
            nodes2 = imat(x-1,:,z); % grab all rows (nodes) of column y-1
            wts = ((1-mat(x,:,z)) + (1-mat(x-1,:,z))) ./ 2;  % average of 1-value for each node
            g = addedge(g,nodes1,nodes2,wts); % add to graph
        end
        if x ~= ydim
            nodes1 = imat(x,:,z);   % grab all rows (nodes) of column y 
            nodes2 = imat(x+1,:,z); % grab all rows (nodes) of column y-1
            wts = ((1-mat(x,:,z)) + (1-mat(x+1,:,z))) ./ 2;  % average of 1-value for each node
            g = addedge(g,nodes1,nodes2,wts); % add to graph
        end
    end
end

for z = 1:zdim  % loop through all "deep" connections in each z plane
    % example:
    %   5  7
    %  /  /
    % 1  3

    for x = 1:xdim
        if z ~= 1
            nodes1 = imat(x,:,z);   % grab all rows (nodes) of column y 
            nodes2 = imat(x,:,z-1); % grab all rows (nodes) of column y-1
            wts = ((1-mat(x,:,z)) + (1-mat(x,:,z-1))) ./ 2;  % average of 1-value for each node
            g = addedge(g,nodes1,nodes2,wts); % add to graph
        end
        if z ~= zdim
            nodes1 = imat(x,:,z);   % grab all rows (nodes) of column y 
            nodes2 = imat(x,:,z+1); % grab all rows (nodes) of column y-1
            wts = ((1-mat(x,:,z)) + (1-mat(x,:,z+1))) ./ 2;  % average of 1-value for each node
            g = addedge(g,nodes1,nodes2,wts); % add to graph
        end
    end
end
disp('done.')

figure
spy(adjacency(g))

答案2:将链接的答案的4个连接的情况扩展到第三维(https://stackoverflow.com/a/3283732/2371031):

disp("calculating adjacency matrix...")

X=3; Y=3; Z=2;
mat = reshape([1:(X*Y*Z)],X,Y,Z);
[x, y, z] = size(mat);                                     % Get the matrix size
diagVec1 = repmat(repmat([ones(y-1, 1); 0], x, 1), z, 1);  % Make the first diagonal vector 
                                                           % (for y-planar connections)
diagVec1 = diagVec1(1:end-1);                              % Remove the last value

diagVec2 = repmat([ones(y*(x-1), 1); zeros(x,1)], z, 1);   % Make the second diagonal vector 
                                                           % (for x-planar connections)
diagVec2 = diagVec2(1:end-x);                              % drop the last x values (zeros)

diagVec3 = ones(x*y*(z-1), 1);                             % make the third diagonal vector 
                                                           % (for z-planar connections)

adj = diag(diagVec1, 1)+diag(diagVec2, y)+diag(diagVec3, x*y);   % Add the diagonals to a zero matrix
adj = adj+adj.';                                           % Add the matrix to a transposed copy of
                                                           %   itself to make it symmetric

disp("done.")

figure
spy(adj)

对于大型情况,Matlab报错:

Error using diag
Requested 65536000x65536000 (32000000.0GB) array exceeds maximum array size preference. Creation of
arrays greater than this limit may take a long time and cause MATLAB to become unresponsive.

Error in blah (line X)
adj = diag(diagVec1, 1) + diag(diagVec2, y) + diag(diagVec3, x*y);   % Add the diagonals to a zero matrix

答案3:一个(非常慢的低效率)for循环版本,它使用权重计算创建一个边列表:

% input matrix of values 
X=3;Y=3;Z=2;
A = reshape([1:(X*Y*Z)],X,Y,Z);

% # nodes is equal to numel(A)
% # edges is equal to number of pairs of 4-connected nodes (cityblock)

nodes1 = [];
nodes2 = [];
wts = [];

xdim = size(A,1);
ydim = size(A,2);
zdim = size(A,3);
% create adjacency matrix
k=1;
disp('finding connected nodes for edge list with weights...')
for z1 = 1:zdim
    for x1 = 1:xdim
        for y1 = 1:ydim
            for z2 = 1:zdim
                for x2 = 1:xdim
                    for y2 = 1:ydim
                        node1 = xdim*(x1-1) + y1 + (xdim*ydim)*(z1-1);
                        node2 = xdim*(x2-1) + y2 + (xdim*ydim)*(z2-1);

                        d = sqrt((x2-x1)^2+(y2-y1)^2+(z2-z1)^2);

                        if d==1

                            % add first node index:
                            nodes1(k) = node1;

                            % add second nodename:
                            nodes2(k) = node2;

                            % add edge weight between these nodes:
                            wts(k) = nanmean([A(x1,y1,z1),A(x2,y2,z2)]);

                            k=k+1;
                        end
                    end
                end
            end
        end
    end
end
g = digraph(nodes1, nodes2, wts);
figure
spy(adjacency(g))

disp('done.')

此人实际上可能适用于大箱子,但我没有耐心寻找。

© www.soinside.com 2019 - 2024. All rights reserved.