MATLAB Answers

0

Replacing for loops for faster computations involving 3D matrix and sampling

Asked by Emily Her on 14 Oct 2018 at 2:59
Latest activity Edited by Bruno Luong on 14 Oct 2018 at 6:45

Hi! I'm hoping there's someone out there who would give me an idea on how to increase computation speed for the following code:

for l=1:size(xSample,2) %%anisotropic 3D gaussian distribution of translation
      for m=1:size(ySample,2) %%anisotropic 3D gaussian distribution of translation
          for n=1:size(zSample,2) %%anisotropic 3D gaussian distribution of translation
              shiftedDoseDist=imtranslate(doseDist,[xSample(l),ySample(m),zSample(n)]); %% 3D image translation
              shiftedDoseDist=imgaussfilt3(shiftedDoseDist,randomError); %% 3D image blurring with Gauss filter
              dose=shiftedDoseDist(index);
              TCP_n(n)=prod(exp(-cellDensity.*vol.*exp(-alpha(q)*numOfFractions*dose- alpha(q)*numOfFractions*dose.^2/AlphaOnBeta+log(2)*T/tPot))); %%function I'm calculating
          end
          TCP_m(m)=1/sum(zWeighting)*dot(zWeighting,TCP_n);
      end
      TCP_l(l)=1/sum(yWeighting)*dot(yWeighting,TCP_m);

So what I'm trying to do is to get a bunch of different combinations of x,y,z from their respective Gaussian distributions for translation of a matrix doseDist (250x250x143), then a Gaussian filter is applied. Then with the resulting distribution, I compute something called TCP where the weighting of the Gaussian distribution is applied to get the expected value.

From the profiler, the bottleneck is in the imtranslate and imgaussfilt3 functions which I don't think they can be further optimised - basically left with getting rid of the nested for loops.

Any help would be appreciated!

  2 Comments

what are the sizes of TCP_n,TCP_m &TCP_l??

There are 3 samples [1x3] for each x,y,z direction and there is another outer for loop which I didn't include (alpha(q)) which is samples [1x4]

Sign in to comment.

1 Answer

Answer by Bruno Luong on 14 Oct 2018 at 6:45
Edited by Bruno Luong on 14 Oct 2018 at 6:45

What are xSample,... zSample and idx?

Not sure you need to translate and filter over and over in the for-loop. It seems awfully wasteful to me.

Rather you could do the filtering once on the big image and extract the information (dose) at appropriate subindexing.

  0 Comments

Sign in to comment.