0

I am currently using fmincon for minimizing a log-Likelihood function in respect to a 18*18 matrix. While on smaller problems the algorithm is very fast, it takes about 2h to converge in the current setup - as I am iterating over this minimisation problem, running through the code may take up to 2 weeks.

Is there a matlab-based, free alternative to fmincon that improves speed on such specific problems? (Costly solutions are discussed here, non-matlab solutions here.) Or would I need to call e.g. a python script from matlab?

The function I want to minimize:

function [L] = logL(A, U, Sigma_e, T, lags)
% A - parameters to optimize w.r.t

logL = 0;
for t = 1 : T - lags
   logL(t, 1) = 0.5*(log(det(A * diag(Sigma_e(t,:)) * A' ) ) + ...
                U(t,:) * (A * diag(Sigma_e(t,:))  * A' )^(-1) *  U(t,:)' );

end
L = sum(logL);

and calling it by:

Options = optimset('Algorithm', 'active-set', 'Display', 'off', 'Hessian','bfgs', ...
'DerivativeCheck','on','Diagnostics','off','GradObj','off','LargeScale','off'); 


A = fmincon( @(A0)logL(A0, U, Sigma_e, T, lags), A0 , [], [] , [] , [] , [] , [] , [], Options); 

(I have tried the different fmincon algorithms without much improvement). Note, T is quite large ~3000. A and A0 are a 18*18 matrices, Sigma_e is T*18, U is T*18

Luks
  • 133
  • 11

1 Answers1

1

I'm not aware of any speedy alternative to fminconst but you can vectorize the logL function to speed up the algorithm. Here is a vectorized version:

function [L] = logL(A, U, Sigma_e, T, lags)
    ia = inv(A);
    iat = ia.';
    N = T - lags;
    UU = zeros(N,1); 
    for t = 1: N
        UU (t) = U(t,:) * (iat .* 1./Sigma_e(t,:) * ia) * U(t,:)';
    end

    L =  0.5 *sum( log(det(A) ^ 2  .* prod(Sigma_e(1:N,:),2)) + UU);
end

In some tests in Octave it nearly is 10X faster than the your solution.

Note that if some elements of Sigma_e are equal to zero you need to compute UU as:

UU (t)=U(t,:) * (A * diag(Sigma_e(t,:))  * A' )^(-1) *  U(t,:)';

These relations are used to convert the loop solution to the vectorized one:

det(a * b * c) == det(a) * det(b) * det(c)

det(a) == det(a.')

det(diag(a)) == prod(a)

(a * b * c)^-1 == c^-1 * b^-1 * a^-1

a * diag(b) == a .* b

inv(diag(a)) == diag(1./a)
rahnema1
  • 15,264
  • 3
  • 15
  • 27
  • Please note: this approach works only for matlab versions 2016 and higher. For versions below, the element-wise multiplications of the first equation of UU(t) would need to be rewritten using for example bsxfun. – Luks May 22 '19 at 16:37