Optimize matlab triple for loop

Question

I have a large matrix samples, I want to make another matrix sampleProb which is the same size, but modified as below:

%samples is a 1000*1000 matrix, M = 1000
sampleProb = zeros(1000);
for i = 1:size(samples, 2)
 for j = 1:M
 for k = 1:M
 if (k ~= j)
 sampleProb(j, i) = sampleProb(j, i) + normcdf(samples(j, i) - samples(k, i));
 end
 end
 end
end
sampleProb = sampleProb./M;

I know it's best to avoid for loops when possible in Matlab and operations should be vectorised, but I can't figure out what to do. How can I optimise this for performance?

• 200_success: Please retitle your post based on the task that your code performs, rather than how you wish to rework it. More background about what you are trying to accomplish would also be appreciated.

Answer 1

Introduction and solution code

Since you mentioned making the output matrix sampleProb, I am assuming it to be initialized as all zeros. So, here's the vectorized implementation with the mighty bsxfun -

%// Get the subtractions
sub1 = bsxfun(@minus,samples,permute(samples,[3 2 1]))
%// Calculate the sum of normcdf's in a vectorized fashion
sampleProb = bsxfun(@minus,sum(normcdf(sub1),3),normcdf(sub1(1,:,1)))./M

---

Benchmarking

Benchmarking Code with M = 100 -

M = 100;
samples = rand(M,M);
sampleProb1 = zeros(M,M);
disp('---------------------------------------- With Original Approach')
tic
for i = 1:size(samples, 2)
 for j = 1:M
 for k = 1:M
 if (k ~= j)
 sampleProb1(j, i) = sampleProb1(j, i) + normcdf(samples(j, i) - samples(k, i));
 end
 end
 end
end
sampleProb1 = sampleProb1./M;
toc, clear sampleProb1 i j k
disp('---------------------------------------- With Proposed Approach')
tic
%// Get the subtractions
sub1 = bsxfun(@minus,samples,permute(samples,[3 2 1]));
%// Calculate the sum of normcdf's in a vectorized fashion
sampleProb = bsxfun(@minus,sum(normcdf(sub1),3),normcdf(sub1(1,:,1)))./M;
toc

Runtime results -

---------------------------------------- With Original Approach
Elapsed time is 21.425729 seconds.
---------------------------------------- With Proposed Approach
Elapsed time is 0.046284 seconds.

• texasflood: Wow, that is fantastic, I did not expect that much of an improvement. Thanks!
• texasflood: The output from both methods produces exactly the same result, save for a tiny difference which is just insignificant quantisation noise (max(max(sampleProb1 - sampleProb))=2.2204e-16) , very impressive! And yes, you were right, sampleProb is initialised to zero, I'll add that to the question