Nested integrals spends quiet different time to calculate for different input values

I'm trying to calculate the cumulative distribution of the functional relationship between several random variables and a specific value's difference (or margin). The mathematical expression should be ∫pdf (x) * pdf (y)* φ (x, y) dxdy. φ The m function calculated for margin.
Assume there are two normal distributions. X follows a normal distribution N(100,10^2), and Y follows a normal distribution N(200, 20^2). I need to calculate the probability of [X+Y ≤ 3/input], where input is a constant variable.
The code I wrote is:
P= integral(@(x) arrayfun(@(x) 1/(2*pi)^0.5/sigma1*exp(-(x-mu1)^2/2/sigma1^2)* integral(@(y) arrayfun(@(y) 1/(2*pi)^0.5/sigma2*exp(-(y-mu2)^2/2/sigma2^2)* margin(x, y, input),y),-inf,+inf),x),-inf,+inf);
function [m] = margin(x, y, input)
I = 3/(x+y);
m = I - input;
if m >= 0
m = 1;
else
m = 0;
end
end
I found that on my computer, when the input is 0.01 and 0.011, the calculation time is 1.099 s and 1597 s, with a significant difference.
What is the reason for the significant difference in calculation time? How to reduce integration calculation time?

Antworten (1)

Torsten
Torsten am 30 Apr. 2024
Bearbeitet: Torsten am 30 Apr. 2024
Check whether the computations are reliable because of the big span of integration.
mu1 = 100;
sigma1 = 10;
mu2 = 200;
sigma2 = 20;
input = 0.011;
P = integral2(@(x,y)normpdf(x,mu1,sigma1).*normpdf(y,mu2,sigma2).*margin(x,y,input),-inf,inf,-inf,inf)
P = 8.4627e-11
function [m] = margin(x, y, input)
I = 3./(x+y);
m = I - input;
for i=1:size(x,1)
for j = 1:size(x,2)
if m(i,j) >= 0
m(i,j) = 1;
else
m(i,j) = 0;
end
end
end
end
E.g. the following integral should give 1, but the result is far off:
mu1 = 100;
sigma1 = 10;
mu2 = 200;
sigma2 = 20;
input = 0.011;
P = integral2(@(x,y)normpdf(x,mu1,sigma1).*normpdf(y,mu2,sigma2),-inf,inf,-inf,inf)
P = 1.2481e-10

Kategorien

Mehr zu Mathematics finden Sie in Hilfe-Center und File Exchange

Gefragt:

am 30 Apr. 2024

Bearbeitet:

am 30 Apr. 2024

Community Treasure Hunt

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

Start Hunting!

Translated by