3
\$\begingroup\$

I am trying to optimise the computation time of my code (for the second time after a first optimisation that gives very good results). Currently, this code is very time-consuming depending on the size of the data (i.e. sometimes very big matrices). So if anyone has some ideas to help me to optimize this code, I will be very grateful!

matrice1=0+(20-0).*rand(500,5);
matrice2=datasample(matrice1,100);

[n1,~]=size(matrice1);
[n2,p2]=size(matrice2); 
% Pre-calculate replicated B and the indices to be modified at each iteration
B_rep=repmat(matrice2,[1 1 n2]);
B_idx=bsxfun(@plus,((0:p2-1)*n2+1)',(0:n2-1)*(n2*p2+1));
B=mean(B_rep,1);
out=zeros(n1,n2); % initialize output array
for i=1:n1
    B_rep(B_idx)=repmat(matrice1(i,:)',1,n2);   
    B_meanB=bsxfun(@minus,B_rep,B); % B minus mean values of B
    A_B_meanB=matrice2'-reshape(B,p2,[]); % A minus B_meanB    

    for j=1:n2
        [~,R]=qr(B_meanB(:,:,j),0);   
        out(i,j)=sum((R'\A_B_meanB(:,j)).^2)*(n2-1);    
    end
end
\$\endgroup\$
3
\$\begingroup\$

Your code looks fairly good, but there are a few things I would do differently:

  • It's very confusing that you call B_rep for B, and call B the mean of B_rep. The comment and code here looks very strange. You should call B_mean = mean(B_rep, 1) to stick with the general naming convention in your code.

    B_meanB=bsxfun(@minus,B_rep,B); % B minus mean values of B
    
  • bsxfun performs better than repmat, so instead of B_rep=repmat(matrice2,[1 1 n2]); you can do:

    B_rep = bsxfun(@times, matrice2, ones(1,1,n2)); 
    
  • Instead of ', you should use .' when transposing an array. The first one is the complex conjugated transpose.

  • i and j are bad variable names in Matlab.

  • You should try to use more spaces. It makes the code much easier to read.


I don't have the statistical toolbox, so I can't test your code myself. I'm not sure how it can be vectorized, so I can't help with much when it comes performance gain I'm afraid.

\$\endgroup\$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.