Design feature matrix of 8-neighborhood elements

Question

I have a matrix (image) of size m x n that I am trying to rearrange differently. I want to design a feature matrix of size (m x n) x 9 where each row corresponds to a matrix element centered around its 8-neighborhood elements. My attempt was to sequentially loop through each element of the original matrix to extract the neighborhood values, however this method is too computationally heavy and takes too long to perform as the matrix size is exhaustively large. Is there anyway to do this cost-beneficially?

Attempt

M_feat = nan(size(img,1)*size(img,2), 9);
temp = zeros(size(img)+2);
temp(2:end-1,2:end-1) = double(img);

for i = 2:size(img,1)+1

    for j = 2:size(img,2)+1

         neighbors = temp(i-1:i+1, j-1:j+1); % read 3-by-3 mini-matrix
         neighbors = reshape(neighbors, 1, []);
         M_feat((i-2)*size(img,1) + (j-1),:) = neighbors; % store row-wise

    end

end

Answer 1

Got it!

new_img = zeros(size(img)+2);
new_img(2:end-1,2:end-1) = double(img);

% Image boundary coordinates without first/last row/column
inner_coord = [2 2; size(new_img,1)-1 size(new_img,2)-1];

% 9x2 matrix with 1 row for the relative shift to reach neighbors
[d2, d1] = meshgrid(-1:1, -1:1);
d = [d1(:) d2(:)];

% Cell array to store all 9 shifted images
temp = {};

for i = 1:size(d,1)
    % x-indices of the submatrix when shifted by d(i,1)
    coord_x = (inner_coord(1,1):inner_coord(2,1)) + d(i,1);
    % y-indices of the submatrix when shifted by d(i,2)
    coord_y = (inner_coord(1,2):inner_coord(2,2)) + d(i,2);
    % image matrix resulting from shift by d(i,)
    temp{i} = reshape(new_img(coord_x, coord_y), 1, []);
end

% Column-wise bind all 9 shifted images (as vectors) from the array
M_feat = vertcat(temp{:}).';