At least one END missing, yet I have written end

12 次查看(过去 30 天)
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!
  7 个评论
Natasha Davina
Natasha Davina 2022-6-8
@Rik okay i will try thank you
@Torsten would this be correct?
y = fsolve(@sumbgg,p,option);

请先登录,再进行评论。

采纳的回答

Jan
Jan 2022-6-8
编辑:Jan 2022-6-8
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.
  3 个评论
Natasha Davina
Natasha Davina 2022-6-8
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
Jan 2022-6-9
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

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Get Started with MATLAB 的更多信息

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by