What's the fastest way to program a matrix recursion in matlab?

The matrix-recursion of the n by n matrices Y_t looks like this:

<<chart-googleapis-com-chart-cht-tx-chl-Y_-t--A-20-2B-20--sum_-i-1---p--20B_-i-Y_-t-i.>>

A and B are given.

This is my attempt, but it runs slowly:

Y = zeros(n,n,T); %Going to fill the 3rd dimension for Y_t, t=1:T
Y(:,:,1:p) = initializingY

for t=(p+1):T
    Y(:,:,t) = A;
    for i=1:p
        Y(:,:,t) = Y(:,:,t) + B(:,:,i)*Y(:,:,t-i);
    end
end

Can you think of a more efficient way to do this?

Answer 1

[Mathias]

Use matrix product!

This code is about 10 times faster on my computer

for t=(p+1):T
    Y(:,:,t) = A + sum(Y(:,:,(t-p+1):t).*B, 3);
end

Regards

[Rosi Marungu]

Unfortunately B also has an index. It is different for the lags from 1 to p. There is no functionality I know of to do matrix multiplication along the 3rd dimension.

[Walter Roberson]

pagefun() can be used to do matrix multiplication along the third dimension