# Multidimensional operations without for loop

2 views (last 30 days)
Aram Eskandari on 27 Mar 2020
Commented: Aram Eskandari on 27 Mar 2020
I have a vector and an input
y = 8.3; % this is different for each run, but not relevant for this question
maxloc = [1 2 1 3 3]; % n long, in this case n = 5, each value goes from 1 to d, in this case d = 3
And a matrix which is attached in .mat-format, it is a double matrix with shape given below
size(theta) = (3,8,5); % d=3, n=5, the columns however are always 8
Which I operate on in the following way
for i=1:n % n = 5 in this case
theta(maxloc(i),7,i) = theta(maxloc(i),7,i) + 0.5;
theta(maxloc(i),8,i) = theta(maxloc(i),8,i) + ((y-theta(maxloc(i),3,i)).^2)*0.5;
theta(maxloc(i),6,i) = 1./(gamrnd(theta(maxloc(i),7,i),1./theta(maxloc(i),8,i)));
theta(maxloc(i),5,i)= theta(maxloc(i),5,i) + 1./theta(maxloc(i),6,i);
theta(maxloc(i),4,i)= (theta(maxloc(i),5,i).*theta(maxloc(i),4,i)+(y./theta(maxloc(i),6,i)) )./(1./theta(maxloc(i),6,i)+theta(maxloc(i),5,i));
theta(maxloc(i),3,i) = normrnd(theta(maxloc(i),4,i),sqrt(1./(theta(maxloc(i),5,i)+1./theta(maxloc(i),6,i))));
end
As you can imagine, the code gets really slow as n increases, so I want to get rid of the for loop.
I haven't included columns 1,2 here because I have a solution for those columns without the use of for loops. However I'm not quite sure how I can get these operations to work in a way which excludes the loops, any help would be greatly appreciated.
EDIT:
The reason this takes so long is because I have 50.000 data-points, meaning y is really 50.000 long and I want to run through my algorithm for each y, and preferably for a large n.

Matt J on 27 Mar 2020
Edited: Matt J on 27 Mar 2020
Q=nan(8,n);
[J,K]=ndgrid(1:8,1:n);
I=maxloc(K);
thetaIndices=sub2ind(size(theta), I,J,K);
QIndices=sub2ind(size(Q), J,K);
Q(QIndices)=theta(thetaIndices);
Q(7,:) = Q(7,:) + 0.5;
Q(8,:) = Q(8,:) + ((y-Q(3,:)).^2)*0.5;
Q(6,:) = 1./(gamrnd(Q(7,:),1./Q(8,:)));
tmp=1./Q(6,:);
Q(5,:)= Q(5,:) + tmp;
tmp2=tmp+Q(5,:);
Q(4,:)= (Q(5,:).*Q(4,:)+(y./Q(6,:)) )./(tmp2);
Q(3,:) = normrnd(Q(4,:),sqrt(1./(tmp2)));
theta(thetaIndices)=Q(QIndices);

Show 1 older comment
Aram Eskandari on 27 Mar 2020
I have edited my question to include the reason why it takes so long to run.
Matt J on 27 Mar 2020
See my edited version, which is loop-less.
Aram Eskandari on 27 Mar 2020
Thanks so much! This seems to be a lot quicker without the for loops, especially with a lot of datapoints. Giving you the green tick.