optimizing nested for loop in matlab

Question

I'm trying to optimize the performance (e.g. speed) of my code. I 'm new to vectorization and tried myself to vectorize, but unsucessful ( also try bxsfun, parfor, some kind of vectorization, etc ). Can anyone help me optimize this code, and a short description of how to do this?

% for simplify, create dummy data
Z = rand(250,1)
z1 = rand(100,100)
z2 = rand(100,100)

%update missing param on the last updated, thanks @Bas Swinckels and @Daniel R
j = 2;
n = length(Z);
h = 0.4;


tic
[K1, K2] = size(z1);
result = zeros(K1,K2);

for l = 1 : K1
    for m = 1: K2
        result(l,m) = sum(K_h(h, z1(l,m), Z(j+1:n)).*K_h(h, z2(l,m), Z(1:n-j)));    
    end
end

result = result ./ (n-j);
toc

The K_h.m function is the boundary kernel and defined as (x is scalar and y can be vector)

function res = K_h(h, x,y)
 res = 0;

 if ( x >= 0 & x < h)
    denominator = integral(@kernelFunc,-x./h,1);  
    res = 1./h.*kernelFunc((x-y)/h)/denominator;
 elseif (x>=h & x <= 1-h)
    res = 1./h*kernelFunc((x-y)/h);
 elseif (x > 1 - h & x <= 1)
    denominator = integral(@kernelFunc,-1,(1-x)./h);
    res = 1./h.*kernelFunc((x-y)/h)/denominator;
 else    
    fprintf('x is out of [0,1]');
    return;
 end
end

It takes a long time to obtain the results: \Elapsed time is 13.616413 seconds.

Thank you. Any comments are welcome. P/S: Sorry for my lack of English

Answer 1

Some observations: it seems that Z(j+1:n)) and Z(1:n-j) are constant inside the loop, so do the indexing operation before the loop. Next, it seems that the loop is really simple, every result(l, m) depends on z1(l, m) and z2(l, m). This is an ideal case for the use of arrayfun. A solution might look something like this (untested):

tic

% do constant stuff outside of the loop
Zhigh = Z(j+1:n);
Zlow = Z(1:n-j);

result = arrayfun(@(zz1, zz2) sum(K_h(h, zz1, Zhigh).*K_h(h, zz2, Zlow)), z1, z2)

result = result ./ (n-j);
toc

I am not sure if this will be a lot faster, since I guess the running time will not be dominated by the for-loops, but by all the work done inside the K_h function.