How to solve the common problem about the using fzero() to solve a integration expression?

2 次查看（过去 30 天）
Chao-Zhen Liu
Chao-Zhen Liu 2018-10-11
评论： Chao-Zhen Liu 2018-10-13
Hi everyone,
I have encountered the problem and get this feedback "Warning: Infinite or Not-a-Number value encountered. "from MATLAB when I use integral() and fzero()...
About this problem, I have tried many advice from my past post on this website but I still could not solve it, hope you can give me some advice.
My code:
tic
U = 14; L = 7; d = (U-L)/2;
T = ( 2*U+L )/(2+1);
Du = U-T; Dl = T - L; d_star = min(Du,Dl);
du = d/Du; dl = d/Dl;
alpha = 0.10;
aLe = 0.05;
N = 2000; % simmulation times
k = 1;
for ksi_given = [-2,-1,-0.5 ,0, 0.5,1,2]
if ksi_given <= 0
sigma = fzero(@(x) aLe - (x^2)*( (dl^2*ksi_given^2+1)/d_star^2 ), 0.4 ) ;
mu = ksi_given*sigma + T;
else
sigma = fzero(@(x) aLe - (x^2)*( (du^2*ksi_given^2+1)/d_star^2 ), 0.4 ) ;
mu = ksi_given*sigma + T;
end
every_sigma(k) = sigma;
every_mu(k) = mu;
j = 1;
for n = [25,50,100,150,200]
i = 1;
for m = 1:1:N
rdata = normrnd(mu,sigma,1,n);
xbar = mean(rdata);
sd = std(rdata);
aLehat_rdata = ( max( (xbar-T)*d/Du,(T-xbar)*d/Dl )/d_star )^2 + var(rdata)/d_star^2;
kursi_max = 0.5; % method1
kursi_hat = (xbar-T)/sd; %method2
fun_max = @(y,C) chi2cdf( (n*( (du.*kursi_max).^2+1 )./C).*aLehat_rdata-y, n-1).* ( 1./ (2.*sqrt(y)) ).* ( (dl^-1).*normpdf( (dl^-1).*sqrt(y) + sqrt(n)*(kursi_max) ) + (du^-1).*normpdf( (du^-1).*sqrt(y) - sqrt(n)*(kursi_max) ) );
exact_UCB_max(i) = fzero(@(C) alpha - integral(@(y) fun_max(y,C),0.00001,(n.*( (du.*kursi_max).^2+1 )./C).*aLehat_rdata), 0.06 );
if kursi_hat <= 0
fun2 = @(y,C) chi2cdf( (n*( (dl.*kursi_hat).^2+1 )./C).*aLehat_rdata-y, n-1).* ( 1./ (2.*sqrt(y)) ).* ( (dl^-1).*normpdf( (dl^-1).*sqrt(y) + sqrt(n)*(kursi_hat) ) + (du^-1).*normpdf( (du^-1).*sqrt(y) - sqrt(n)*(kursi_hat) ) );
% plotfun = @(C) alpha - integral(@(y) fun2(y,C),0.00001,(n.*( (dl.*kursi_hat).^2+1 )./C).*aLehat_rdata);
% c = linspace(1e-2,0.2,100);
% for i=1:numel(c)
% fc(i) = plotfun(c(i));
% end
% plot(c,fc)
exact_UCB_hat(i) = fzero(@(C) alpha - integral(@(y) fun2(y,C),0.00001,(n.*( (dl.*kursi_hat).^2+1 )./C).*aLehat_rdata), 0.06 );
else
fun2 = @(y,C) chi2cdf( (n*( (du.*kursi_hat).^2+1 )./C).*aLehat_rdata-y, n-1).* ( 1./ (2.*sqrt(y)) ).* ( (dl^-1).*normpdf( (dl^-1).*sqrt(y) + sqrt(n)*(kursi_hat) ) + (du^-1).*normpdf( (du^-1).*sqrt(y) - sqrt(n)*(kursi_hat) ) );
exact_UCB_hat(i) = fzero(@(C) alpha - integral(@(y) fun2(y,C),0.00001,(n.*( (du.*kursi_hat).^2+1 )./C).*aLehat_rdata), 0.06 );
end
each_aLehat(i) = aLehat_rdata;
i = i + 1;
end
CR_1(j,k) = mean(aLe<exact_UCB_max); %method1
CR_2(j,k) = mean(aLe<exact_UCB_hat); %method2
average_aLehat(j,k) = mean(each_aLehat);
sd_aLehat(j,k) = std(each_aLehat);
average_exact_UCB_max(j,k) = mean(exact_UCB_max);
sd_exact_UCB_max(j,k) = std(exact_UCB_max);
average_exact_UB_hat(j,k) = mean(exact_UCB_hat);
sd_exact_UB_hat(j,k) = std(exact_UCB_hat);
j = j + 1;
end
k = k + 1;
end
toc
Error message:
Warning: Infinite or Not-a-Number value encountered.
> In integralCalc/iterateScalarValued (line 349)
In integralCalc/vadapt (line 132)
In integralCalc (line 75)
In integral (line 88)
In @(C)alpha-integral(@(y)fun_max(y,C),0.00001,(n.*((du.*kursi_max).^2+1)./C).*aLehat_rdata)
In fzero (line 522)
About this problem, I am very sure that the problem will be caused by the two lines code below,
exact_UCB_max(i) = fzero(@(C) alpha - integral(@(y) fun_max(y,C),0.00001,(n.*( (du.*kursi_max).^2+1 )./C).*aLehat_rdata), 0.06 );
exact_UCB_max(i) = fzero(@(C) alpha - integral(@(y) fun_max(y,C),0.00001,(n.*( (du.*kursi_max).^2+1 )./C).*aLehat_rdata), 0.06 );
because I have encountered the similar problem before, so I try to modify the lower bound of integral from 0 to 0.00001 and use the plot technique like below to decide the starting value of fzero().
plotfun = @(C) alpha - integral(@(y) fun2(y,C),0.00001,(n.*( (dl.*kursi_hat).^2+1 )./C).*aLehat_rdata);
c = linspace(1e-2,0.2,100);
for i=1:numel(c)
fc(i) = plotfun(c(i));
end
plot(c,fc)
Finally, I use 0.06 as the starting value for fzero() but it still return the same error message ...
Does anyone have any idea of this problem to let me fix it?
Very thank you!
  7 个评论
显示 5更早的评论
Walter Roberson
Walter Roberson 2018-10-12
Can you tell me how to evaluate the fun_max value like what you do?
It involves my adding symbolic replacements for the functions such as chi2cdf that you use, so that I can call the functions with symbolic numbers and get a more precise result than the nan or inf that you are getting returned. Unfortunately there are parts of the existing numeric statistical functions that do not accept symbolic inputs, because they try to do things like
a(a < 0) = nan;
to try to deal with values out of range. This is a problem if the inputs involve symbolic variables, because symbolic_variable relation constant is not decideable without a particular value for the variable -- so I have to rewrite in terms of detecting symbolic variables and using piecewise() to do the tests.
... it is a bit slow going to make the changes.
Chao-Zhen Liu
Chao-Zhen Liu 2018-10-13
Hi Walter,
Thank for your help and explanation!
But it seems like my integration expression could not work by numerical integration and it is also have to do many pretreatment if we want to use symbolic integration, right?

请先登录，再进行评论。

请先登录，再回答此问题。

回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Calculus 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by