Problem finding minimum of a function

Hello,

I have a function f = (x,n,mu1,mu2,sigma1,sigma2,prop2), which is a mixed Gumbel distribution. I want to estimate the parameters mu1,mu2,sigma1,sigma2,prop2 dependend on my collected data x and n. I'm trying to find the values of the parameters which minimize L = sum(log(f = (x,n,mu1,mu2,sigma1,sigma2,prop2))) using fminunc or fminsearch.

When i try to run my code there is an error:

Not enough input arguments.

Error in max_like>myfunc (line 25)

L = sum(log(f(x_dat,n_dat,mu1,mu2,sigma1,sigma2,prop2)));

Error in fminsearch (line 200)

fv(:,1) = funfcn(x,varargin{:});

Error in max_like (line 10)

x = fminsearch(@myfunc,[25, 35, 5, 5, 0.5]);

When I only call the function myfunc it works fine. It has something to do with evaluating the function f() inside myfunc(). I also tried using mle(), but this doesnt work with two variables. I know using global variables isn't a good practise, but this function will be implemented into a GUI, where I have private properties for the whole program and I don't need global variables there.

Thanks for your help in advance.

global x_dat;
x_dat = [12;13;14.5;15];
global n_dat;
n_dat = [2;2;5;5];

%x = fminunc(@myfunc,[16, 12, 5, 5, 0.5]);
x = fminsearch(@myfunc,[16, 12, 5, 5, 0.5]);

function L = myfunc(mu1,mu2,sigma1,sigma2,prop2)

    global x_dat
    global n_dat
    F1 = @(x,mu1,sigma1)(exp(-exp(-(x-mu1)/sigma1)));
    F2 = @(x,mu2,sigma2)(exp(-exp(-(x-mu2)/sigma2)));
    f1 = @(x,mu1,sigma1) ((1/sigma1)*exp(-((x-mu1)/sigma1 + exp(-(x-mu1)/sigma1)))).* (exp(-exp(-(x-mu1)/sigma1)));
    f2 = @(x,mu2,sigma2) ((1/sigma2)*exp(-((x-mu2)/sigma2 + exp(-(x-mu2)/sigma2)))).* (exp(-exp(-(x-mu2)/sigma2)));
    f_mix = @(x,mu1,mu2,sigma1,sigma2,prop2)(1-prop2)*f1(x,mu1,sigma1)+ prop2*f2(x,mu2,sigma2);
    
    f = @(x,n,mu1,mu2,sigma1,sigma2,prop2)(n.*((1-prop2).*F1(x,mu1,sigma1) + prop2.*F2(x,mu2,sigma2)).^(n-1)).*f_mix(x,mu1,mu2,sigma1,sigma2,prop2);
    
    L = sum(log(f(x_dat,n_dat,mu1,mu2,sigma1,sigma2,prop2))); 
       
end

Answer 1

[Stephan]

Hi,

try this:

[x, fval] = curvefit_gumble

function [x,fval] = curvefit_gumble
x_dat = [12;13;14.5;15];
n_dat = [2;2;5;5];
oldopts = optimset('fminsearch');
options = optimset(oldopts,'MaxFunEvals',250,'MaxIter',250,'Display','Off');
[x,fval] = fminsearch(@myfunc,[16, 12, 5, 5, 0.5],options);

function L = myfunc(x)
mu1 = x(1);
mu2 = x(2);
sigma1 = x(3);
sigma2 = x(4);
prop2 = x(5);
F1 = @(x,mu1,sigma1)(exp(-exp(-(x-mu1)/sigma1)));
F2 = @(x,mu2,sigma2)(exp(-exp(-(x-mu2)/sigma2)));
f1 = @(x,mu1,sigma1) ((1/sigma1)*exp(-((x-mu1)/sigma1 + exp(-(x-mu1)/sigma1)))).* (exp(-exp(-(x-mu1)/sigma1)));
f2 = @(x,mu2,sigma2) ((1/sigma2)*exp(-((x-mu2)/sigma2 + exp(-(x-mu2)/sigma2)))).* (exp(-exp(-(x-mu2)/sigma2)));
f_mix = @(x,mu1,mu2,sigma1,sigma2,prop2)(1-prop2)*f1(x,mu1,sigma1)+ prop2*f2(x,mu2,sigma2);
f = @(x,n,mu1,mu2,sigma1,sigma2,prop2)(n.*((1-prop2).*F1(x,mu1,sigma1) + prop2.*F2(x,mu2,sigma2)).^(n-1)).*f_mix(x,mu1,mu2,sigma1,sigma2,prop2);
L = sum(log(f(x_dat,n_dat,mu1,mu2,sigma1,sigma2,prop2)));
end
end

Results:

Exiting: Maximum number of iterations has been exceeded
         - increase MaxIter option.
         Current function value: -47.742019 


x =

   23.3458   13.8903    5.3200    6.4483    0.1364


fval =

  -47.7420

Increasing MaxFunEvals and MaxIter does not change the result.

[snipsnap333]

Thank you very much!