Search All of the Math Forum:
Views expressed in these public forums are not endorsed by
Drexel University or The Math Forum.
|
|
|
|
MLE - use with Discrete Weibull mean function
Posted:
Nov 7, 2009 1:31 PM
|
|
I am having trouble with an MLE function that was passed along to me. I am new to MATLAB and am having trouble shooting problems. I am taking 'v' values that are randomly generated that are discrete numbers that range from 0 to a couple thousand. In the code, "v' can have a record length of 10, 50, 100, 1000, and 10,000. In the Llikef function I would like to take the log of the function, but I have been getting log zero when large numbers usually over 1,022 are entered into the function. So I took the log out. The thing I need help with is that I want my q and N (eta) values to be positive. The function sometimes outputs negative q's and eta's. Does anybody know how to do this? Thank you in advance!
In the function that(1) = q and that(2) = eta.
vsim are whole numbers that range from 0 to a few thousand
%find the Mean MLE Fitted Average.
%1) take each record length of 'v' values and plug it into the MLE routine
%2) use the q and eta values generated from that function and plug them
%in the DW mean function
slength=1;
elength=recordlength;
meanindex=1;
vsim;
while meanindex<=100
if sum(vsim(slength:elength))>0
x=vsim(slength:elength);
Llikef = @(that)-sum(max(that(1),.00001).^(x.^max(that(2),.00001))-max(that(1),.00001).^((x+1).^max(that(2),.00001)));
%Llikef = @(that)-sum(log(max(that(1),.00001).^(x.^max(that(2),.00001))-max(that(1),.00001).^((x+1).^max(that(2),.00001))))
that0 = [0.5 1];
[that,Llikef] = fminsearch(Llikef,that0);
if Llikef<=0
that(1)=.5;
that(2)=1;
else end
if that(1)>0
q2(meanindex)=that(1);
else q2(meanindex)=.5;
end
if that(2)>0;
eta2(meanindex)=that(2);
else eta2(meanindex)=1;
end
N=eta2(meanindex);
q=q2(meanindex);
M=1e5;
vm=1:M;
mmlefm(meanindex)= sum(q.^(vm.^N))+(gamma(1/N)*gammainc((M+1)^N*(-log(q)),1/N,'upper'))/(N*(-log(q))^(1/N));
meanindex=meanindex+1;
slength=elength+1;
elength=elength+recordlength;
else
q2(meanindex)=.5;
eta2(meanindex)=1;
mmlefm(meanindex)=0;
slength=elength+1;
elength=elength+recordlength;
meanindex=meanindex+1;
end
end
|
|
|
|
