3
votes

in the MATLAB code below, the matrix B is created from matrix A by randomly replacing a portion of its elements (empty_x) with 6. Instead of replacing the values of empty_x elements randomly throughout the matrix, how could I keep them clumped (but still random placement of the clump)?

That is I want all of the 6-value elements to be neighboring. If it easiest to replace empty_x elements a the square matrix A with smaller, subset-ed square matrix (size empty_x) to create matrix B then that would be okay. It would be cool, but not necessary, to have the clumps that were not always square matrices (i.e heterogeneity).

I would appreciate some ideas on how to accomplish this task.

Cheers.

A = [1 2 3 4 5; 5 4 3 2 1; 1 4 3 2 5; 4 3 2 1 5; 2 1 3 5 4];
B = A;
nA = numel(A);
empty_x = randi(10);
B(randperm(nA,(empty_x))) = 6;
2
So you've got a square matrix, and you want some random square subset of that matrix to be replaced with 6's? Is that the problem? I'm not sure what you mean by 'clumped'Nato Saichek
@NatoSaichek That would work. I was using clumped to mean that the elements with 6's would touch each other as in a subset square matrix, or some other shape. Thanks.nofunsally
Take a look at this linkbonCodigo
@nofunsally What other constraints are you envisioning? Does every matrix need to have the same number of 6s inserted? Do the clumps need to be uniformly distributed throughout the matrix? Can they touch?jerad
@jerad For any given simulation the number of 6's to be inserted can be the same. Clumps do not need to be uniformly distributed. @bonCodigo I am looking at MATLAB's datasample now. Thanks.nofunsally

2 Answers

1
votes

My approach would be the following:

1) Generate a single random number (uniform distribution)  
     on the interval `[1 numel(A)]`. Use this as the linear index  
     of a seed for your clump.  

while clump_size < desired_clump_size     
    2) Generate a list of all positions in the matrix adjacent to  
         (but not already included in) the existing clump.  
    3) Randomly select one of these indices  
    4) Grow the clump by placing an element in this position.  
end 

I'm not going to write the code; it shouldn't be difficult to implement, particularly if this piece of code isn't a performance bottleneck in your overall project.

EDIT: Since you gave it a try yourself, here's some code:

desired_clump = 5;
matrix_size = 5;
A = zeros(matrix_size);
[C,R]=meshgrid(1:size(A,1), (1:size(A,2))'); %'# row and column numbers for each element

seed = ceil(rand(1)*numel(A));
#% I would have used randi(1) but octave online utility doesn't have it
A(seed) = 1; #% initialize a clump
clump_size = 1;

while clump_size < desired_clump
    CI = A==1; #% logical index of current clump
    CR = reshape(R(CI),1,1,[]); #% 1x1xN index of row values of current clump
    CC = reshape(C(CI),1,1,[]); #% 1x1xN index of col values of current clump
    ADJ = sum(bsxfun(@(x,y)abs(x-y),R,CR)<=1 & bsxfun(@(x,y)abs(x-y),C,CC)<=1, 3)>0 & ~A;
    #% ADJ is the indices of the elements adjacent to the current clump
    B=A; #% for display purposes only
    B(ADJ)=2;
    disp(B)
    disp(' ')
    POS = find(ADJ); #% linear indices of the adjacent elements
    IND = ceil(rand(1)*numel(POS)); #% random index into POS vector
    A(POS(IND))=1; #% grow the clump
    clump_size = clump_size+1;
end
disp(A);

Output: 1 indicates elements in the clump; 2 means eligible for clump expansion

iteration 1:
   0   0   2   1   2
   0   0   2   2   2
   0   0   0   0   0
   0   0   0   0   0
   0   0   0   0   0
iteration 2: 
   0   0   2   1   2
   0   0   2   1   2
   0   0   2   2   2
   0   0   0   0   0
   0   0   0   0   0
iteration 3: 
   0   0   2   1   2
   0   0   2   1   2
   0   0   2   1   2
   0   0   2   2   2
   0   0   0   0   0
iteration 4: 
   0   0   2   1   1
   0   0   2   1   2
   0   0   2   1   2
   0   0   2   2   2
   0   0   0   0   0

Final clump:
   0   0   0   1   1
   0   0   1   1   0
   0   0   0   1   0
   0   0   0   0   0
   0   0   0   0   0

Generating a single random number each time shouldn't be all that slow. There are also undoubtedly ways to speed it up if it is truly a bottleneck. Hopefully this example can get you a little further though.

0
votes

Using some of the tips above I constructed the code below. It works fine for smaller sizes of A, but it is incredibly slow when A is large. Thanks for the guidance.

clear all
A = zeros(40,40);
[M N] = size(A);
B = A;
nA = numel(A);
per_clump = 10;
dClump = nA*(per_clump/100);
seed = randi(nA);
clumpers = zeros(8,1);
new_seed = seed;
counter = 0;
while counter < dClump; % size of clump
    seed = new_seed;

    for iSize = seed; % find adjacent elements

        west = iSize - M;
        if west < 1
            west = iSize; % the boundary is not periodic
        end;

        east = iSize + M;
        if east > nA
            east = iSize;
        end; %

        north = iSize - 1;
        if north < 1
            north = iSize;
        end; %

        south = iSize + 1;
        if south > nA
            south = iSize;
        end; %

        nwest = iSize - M - 1;
        if nwest < 1
            nwest = iSize;
        end; %

        neast = iSize + M - 1;
        if neast > nA
            neast = iSize;
        end; %

        swest = iSize - M + 1;
        if swest < 1
            swest = iSize;
        end; %

        seast = iSize + M + 1;
        if seast > nA
            seast = iSize;
        end; %

        clumpers = [(west) (east) (north) (south) (nwest) (neast) (seast) (swest)]; % index of adjacent elements
        %new_seed = randsample(clumpers,1); % pick one, really slow

        z = randperm(size(clumpers,2)); % this also really slow
        new_clumpers = clumpers(z);
        new_seed = new_clumpers(randi(8));

        if B(new_seed) == 6;
            %B(new_seed) = B(seed);
            new_seed = seed;
            counter = counter;
        else
            B(new_seed) = 6;
            counter = counter+1;
        end;
    end; %end adj element
end; % end clump size