21 views (last 30 days)

I am kind of new to MATLAB and I want to obtain the empitical cumulative distribution function (CDF) of the below PDF:

where:

the PDF function is as follows:

GEV_function=@(y) (Gamma_d.*Gamma_D)./(Delta_t.*Etta_d.*Etta_D).*(1./y).*(exp(-(1-(Mu_d.*Gamma_d)+(Gamma_d.*x)./y).^((-1)./Etta_d))./exp((1-(Mu_D.*Gamma_D)+(Gamma_D.*y)./Delta_t).^((-1)./Etta_D))).*(1-(Mu_d.*Gamma_d)+(Gamma_d.*x)./y).^(-1-(1./Etta_d)).*(1-Mu_D.*Gamma_D+Gamma_D.*y./Delta_t).^(-1-(1./Etta_D));

PDF=integral(GEV_function,0,inf);

To obtain the CDF and to check if I get CDF = 1 as x approuches to infinity, I have integrated the f(x) from minus infinity to positive infinity as follows:

GEV_Original=@(y,x) (Gamma_d.*Gamma_D)./(Delta_t.*Etta_d.*Etta_D).*(1./y).*(exp(-(1-(Mu_d.*Gamma_d)+(Gamma_d.*x)./y).^((-1)./Etta_d))./exp((1-(Mu_D.*Gamma_D)+(Gamma_D.*y)./Delta_t).^((-1)./Etta_D))).*(1-(Mu_d.*Gamma_d)+(Gamma_d.*x)./y).^(-1-(1./Etta_d)).*(1-Mu_D.*Gamma_D+Gamma_D.*y./Delta_t).^(-1-(1./Etta_D));

CDF = integral2(GEV_Original,0,inf,-inf,inf)

But the answer is 0.6449 and not 1.

I guess the usage of double integral is leading to such an error. Is my code correct?

I appreciate any thoughts or help.

David Goodmanson
on 18 Dec 2019

Hello Manesf,

This pdf seems related to the Weibull distribution. I made some abrreviations

Etta_d = xid Etta_D = xiD

Gamma_d = gd Gamma_D = gD

Mu_d = mud Mu_D = muD,

and changed the sign of the argument of the second exp, multiplying by it rather than dividing. This leads to

PDF = @(x,y) (gd*gD)/(Del*xid*xiD)*(1./y) ...

.*exp(-(1 -mud*gd +gd*x./y ).^(-1/xid)).*(1 -mud*gd +gd*x./y ).^(-1-1/xid) ... [1]

.*exp(-(1 -muD*gD +gD*y/Del).^(-1/xiD)).*(1 -muD*gD +gD*y/Del).^(-1-1/xiD); [2]

Take a look at line [2] , considered as a distribution on its own. There is a basic normalized distribution (whose name I don't know),

(1/xiD) * exp(-y.^(-1/xiD)) .* y.^(-1-1/xiD)

with 0 <= y <= inf. Variable y cannot be negative since it is taken to a fractional power.

Line [2] contains a scaled and shifted version of y. Now the argument itself must be positive,

0 <= (1 -muD*gD +gD*y/Del) <= inf

The net result is that y can take on some negative values, and the mean value of y turns out to be y_mean = 1, which I assume is intentional. The lower limit for y is

y0 = ((muD*gD-1)*Del/gD)

which is negative. If you integrate line [2] from y0 to inf, and include the normalization factors in front that are associated with that line,

gD/(Del*xiD)

you get 1.

For the full pdf, integration over y is complcated because now you have the condition in line [1]

0 <= (1 -mud*gd +gd*x./y ) <= inf

and the y limits depend on x. But you sensibly wanted to do both the x integration and the y integration as a check, which is fairly easy. Substitute

x = y*z dx = y*dz

and do the z integration first. The y in y*dz cancels the 1/y in the pdf and the resulting integration is from z0 to inf, where z0 can be calculated in the same way that y0 was. The result is

Del = 25;

xid = .02;

xiD = .08;

mud = -2.95;

muD = .05;

gD = (gamma(1-xiD)-1)/(1-muD);

gd = (gamma(1-xid)-1)/(1-mud);

PDF1 = @(z,y) (gd*gD)/(Del*xid*xiD) ...

.*exp(-(1 -mud*gd +gd*z ).^(-1/xid)).*(1 -mud*gd +gd*z ).^(-1-1/xid) ...

.*exp(-(1 -muD*gD +gD*y/Del).^(-1/xiD)).*(1 -muD*gD +gD*y/Del).^(-1-1/xiD);

z0 = ((mud*gd-1)/gd)*.9999;

y0 = ((muD*gD-1)*Del/gD)*.9999;

CDF1 = integral2(PDF1,z0,inf,y0,inf)

format long

CDF1 =

0.999999898129661

There are some integration issues near z0 and y0. To address that, I just multiplied each limit by .9999. To really do it right would take some work but I think the simple approach here shows that the result is correct.

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

Start Hunting!
## 9 Comments

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778274

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778274

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778279

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778279

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778315

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778315

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778540

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778540

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778543

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778543

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778587

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778587

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778589

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778589

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778592

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778592

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778619

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/496802-how-to-obtain-cdf-from-the-below-pdf-function#comment_778619

Sign in to comment.