Vectorize a doble loop for

Hi everyboy

I am trying to vectorize the following code, but I can't. I want to calculate the matrix p_hat from the matrixs p and omega (these last matrix are of the same order) . The vector a is a function that gives me another index.

Could anyone help me, please?

for k=1:length(a)
    for j=1:length(a)
        if omega(j,k)==1    
            p_hat(k,a(j))=p(j,k);
        end
    end
end

[Stephen23]

" The vector a is a function that gives me another index"

It can't be both: is it a vector, or a function?

[jose ramirez alvarez]

It is a vector

Answer 1

[Jan]

Assuming that omega is a logical matrix (please post what the inputs are, such that the readers do not have to guess the details):

p_hat = zeros(length(a), length(a));  % Needed only if p_hat is existing before
p_hat(:, a) = (p .* omega).';

Some assumed test data:

omega = rand(10, 10) > 0.2;
a     = randperm(10);
p     = rand(10, 10);

p_hat2(:, a) = (p .* omega).';        % Vectorized

p_hat = zeros(length(a), length(a));  % Loop
for k=1:length(a)
   for j=1:length(a)
      if omega(j,k)==1    
          p_hat(k,a(j))=p(j,k);
      end
   end
end

isequal(p_hat, p_hat2)                 % Compare
>> logical(1)

[jose ramirez alvarez]

Yes, it works. But how to vectorize de outer loop? Maybe something like

    p_hat(1:length(a), a_aux(omega(:, 1:length(a)))) = p(omega(:, 1:length(a)), 1:length(a));

[Jan]

See the edited answer.

[jose ramirez alvarez]

Thanks a lot!

Answer 2

[Bruno Luong]

p = rand(3,4)
a = ceil(10*rand(size(p,1),1))
omega = rand(size(p))>0.4

[j,k] = find(omega)
p_hat = accumarray([k,a(j)],p(omega))