At least one END missing, yet I have written end

Hi, I am trying to estimate the parameters of a distribution using the maximum likelihood estimator, and this is my loglikelihood function/objective function

I have written a code for the formula for the objective function as so (filename sumbgg)

function F = sumbgg(p,x)

sum1 = 0;

sum2 = 0;

sum3 = 0;

sum4 = 0;

sum5 = 0;

sum6 = 0;

sum7 = 0;

sum8 = 0;

sum9 = 0;

sum10 = 0;

for i=1:length(x)

sum1 = sum1 + log(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)));

sum2 = sum2 + ((1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))^p(1))*log(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))/log(1-(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))^p(1));

sum3 = sum3 + exp(p(3)*x(i));

sum4 = sum4 + (exp(p(3)*x(i)-1)*exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))/(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)));

sum5 = sum5 + ((exp(p(3)*x(i)-1)*exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))*(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))^(p(1)-1))/(1-(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))^p(1));

sum6 = sum6 + x(i);

sum7 = sum7 + (exp(p(3)*x(i)))*(1-p(3)*x(i));

sum8 = sum8 + (exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))*(p(3)*x(i)*exp(p(3*x(i)))+exp(p(3)*x(i))+1)/(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)));

sum9 = sum9 + (((1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))^p(1))*exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))*(p(3)*x(i)*exp(p(3*x(i)))+exp(p(3)*x(i))+1)/(1-(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))^p(1));

sum10 = sum10 + log(1-((1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1))))^p(1));

end

F = [(length(x)/p(1))+(p(4)*sum1)-((p(5)-1)*sum2);

(length(x)/p(2))+(length(x)/p(3))-(1/p(3))*sum3+(((p(1)*p(4)-1)/p(3))*sum4)+((p(1)*(p(5)-1)/p(3))*sum5);

(-length(x)*p(2)/(p(3)^2))+ sum6 + p(2)*sum7/(p(3)^2)-p(2)*(p(4)*p(1)-1)*sum8/(p(3)^2)+p(1)*p(2)*(p(5)-1)*sum9/(p(3)^2);

-length(x)*(psi(p(4))+psi(p(4)+p(5)))+p(1)*sum1;

-length(x)*(psi(p(4))+psi(p(4)+p(5)))+p(1)*sum10];

F = double(F);

end

And I want to find the values for each parameter as below (filename calc)

p = [2 0.005 1.16 0.01 0.01];

x = xlsread('Time of First Failure of 500 MW Generators.xlsx');

i = 1;

while(max(sumbgg(p,x))>=10^4 && i<=9999)

option = optimset();

y = fsolve(@sumbgg,p,option,x);

p = y;

fprintf('Iterasi ke-%d\n',i);

i = i + 1;

end

y = sumbgg(x,a);

When I tried to run it returned this error:

Error: File: calc.m Line: 5 Column: 1

At least one END is missing. The statement beginning here does not have a matching end.

I don't understand what I have to do since I have ended the while loop, does anybody know how to fix this?

Also, I have read about optimset but I still don't understand what method it uses to optimize the functions is it perhaps Newton Raphson?

I am new to Matlab so any advice would be much appreciated, thank you!

[Natasha Davina]

slight edit the last line of the code shouldve been

y = sumbgg(p,x);

[Rik]

I suspect you have mismatched parentheses. I would suggest using more whitespace and line continuations (...) to make it more readable. That should also help to make sure everything is matched as it should be.

[Stephen23]

@Natasha Davina: please upload the Mfile by clicking the paperclip button.

[Torsten]

@Natasha Davina

Most probably not the reason for the error message, but nevertheless an error: Your call to "fsolve".

[Jan]

I get a different error message:

Array indices must be positive integers or logical values.
Error in bbsd (line 21)
    sum8 = sum8 +
    (exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)))*(p(3)*x(i)*exp(p(3*x(i)))+exp(p(3)*x(i))+1)/(1-exp(-(p(2)/p(3))*(exp(p(3)*x(i))-1)));     

See my answer for the fix of assumed typo.

Error: File: calc.m Line: 5 Column: 1

Please tell us, which one the line 5 is. Is it:

    option = optimset();    

If so, is this the original optimset function? By the way, move this before the loop to save processing time.

[Natasha Davina]

@Stephen23 here are the files i used

[Natasha Davina]

@Rik okay i will try thank you

@Torsten would this be correct?

y = fsolve(@sumbgg,p,option);

Answer 1

[Jan]

At least one bug is here:

sum8 = sum8 + (exp(-(p(2) / p(3)) * (exp(p(3) * x(i)) - 1))) * ...
    (p(3) * x(i) * exp(p(3 * x(i))) + exp(p(3) * x(i)) + 1) / ...
    ... %                ^^^^^^^^^  This is no valid index for p
    ... %          exp(p(3) * x(i)) is meant
    (1 - exp(-(p(2) / p(3)) * (exp(p(3) * x(i)) - 1)));

This appears again:

sum9 = sum9 + (((1 - exp(-(p(2) / p(3)) * ...
    (exp(p(3) * x(i)) - 1)))^p(1)) * exp(-(p(2) / p(3)) * ...
    (exp(p(3) * x(i)) - 1))) *  ...
    (p(3) * x(i) * exp(p(3*x(i))) + exp(p(3) * x(i))+1) / ...
    ... %               ^^^^^^^^
    (1 - (1 - exp(-(p(2)/p(3)) * (exp(p(3) * x(i)) - 1))) ^ p(1));

I think, this is much easier to debug:

function F = sumbgg(p,x)
sum1 = 0;
sum2 = 0;
sum3 = 0;
sum4 = 0;
sum5 = 0;
sum6 = 0;
sum7 = 0;
sum8 = 0;
sum9 = 0;
sum10 = 0;
n = length(x);
for i = 1:n
   e = exp(p(3) * x(i));
   d = exp(-p(2) / p(3) * (e - 1));
   c = 1 - d;
   b = exp(p(3) * x(i) - 1);
   a = c^p(1);
   
   sum1 = sum1 + log(c);
   sum2 = sum2 + a * log(c) / log(1 - a);
   sum3 = sum3 + e;
   sum4 = sum4 + b * d / c;
   sum5 = sum5 + b * d / c * a / (1 - a);
   sum6 = sum6 + x(i);
   sum7 = sum7 + e * (1 - p(3)*x(i));
   sum8 = sum8 + d * (p(3)*x(i) * e + e + 1) / c;
   sum9 = sum9 + (a * d) * (p(3) * x(i) * e + e + 1) / (1 - a);
   sum10 = sum10 + log(1 - a);
end

n = length(x);
p3_2 = p(3)^2;
F = [n / p(1) + p(4) * sum1 - (p(5) - 1) * sum2;
     n / p(2) + n / p(3) - 1 / p(3) * sum3 + ...
     (p(1) * p(4) - 1) / p(3) * sum4 + p(1) * (p(5) - 1) / p(3) * sum5;
     -n * p(2) / p3_2 + sum6 + p(2) * sum7 / p3_2 - ...
     p(2) * (p(4) * p(1) - 1) * sum8 / p3_2 + ...
     p(1) * p(2) * (p(5) - 1) * sum9 / p3_2;
     -n * (psi(p(4)) + psi(p(4) + p(5))) + p(1) * sum1;
     -n * (psi(p(4)) + psi(p(4) + p(5))) + p(1) * sum10];
F = double(F);
end

Remove clutter from the code: Unnecessary parentheses (if they do not support clarity), repeated expensive calculations (exp and power are very expensive!). USe temporary variables for repeated calculations. Then the inner structure of the formula gets more obvious and this supports the debugging.

[Natasha Davina]

thank you

[Natasha Davina]

Thank you very much Jan, I ran the code and it was able to run!! If I may ask a followup question, this is the outcome that I got

Can I ask what is 1.0e+03*?? And do you have any idea how to stop getting NaN values? Thanks!

[Jan]

Matlab displays the values of arrays in a scaled format. 1.0e+03 mean 1.0 * 10^3, or 1000.

Usually NaN are produced by dividing by zero. Simply let Matlab stop, if this happens and check the local variables:

dbstop if naninf