Integrating small valued functions

2 views (last 30 days)
Paulo F.
Paulo F. on 17 Feb 2015
Commented: Paulo F. on 20 Feb 2015

Well, maybe it is better to start with the mathematical part of the problem:

I want to obtain the expected highest number among N random numbers generated by a Normal (m,s).

The expression for this would be (in LaTeX notation), f() and F() are density and cumulative functions respectively:

$$\frac{\int_{-\infty}^\infty x f(x|m,s) F^{n-1}(x|m,s)dx}{\int_{-\infty}^\infty f(x|m,s) F^{n-1}(x|m,s)dx}$$

Now, for big moments (m,s) this gives crazy values, and big moments I refer to things like 30 and above. Also crazy things start to happen with values of N bigger than 35.

An example of what happens is giving a result below the mean, when it is easy to proof that this number will always be bigger than the mean (note F()<1, => F()^k = over-representing the numbers to right of the distribution).

In code, what I do is to set up the parameters and then define the function:

f1=@(x)x.*pdfnorm(x,m,s).*cdfnorm(x,m,s).^(n-1);
f1s=@(x)pdfnorm(x,m,s).*cdfnorm(x,m,s).^(n-1);
p1=integral(@(x)f1s,-Inf,Inf);
p2=integral(@(x)f1,-Inf,Inf);
result=p2/p1;

I have read that Integral might have some troubles with very small or big numbers, as it would be the case for $x$ at some point, or cdfnorm()^{n-1} with high $n$ values. I have not find a way to solve this, do you have any clue?.

I thank you all in advance,

Cheers!

  5 Comments
Paulo F.
Paulo F. on 20 Feb 2015
Thanks for your answers, I will try to do it in the expression instead of using the pdf and cdf functions to see if it works. I will let you know.
Thanks again

Sign in to comment.

Accepted Answer

Mike Hosea
Mike Hosea on 20 Feb 2015
Could you try this and tell me what values of m, s, and n lead to unexpected results? I added some waypoints near the center of the distribution for fear that sometimes the integrator might "miss the action" with its initial mesh.
function result = testfun(m,s,n)
f1=@(x)x.*normpdf(x,m,s).*normcdf(x,m,s).^(n-1);
f1s=@(x)normpdf(x,m,s).*normcdf(x,m,s).^(n-1);
w = m + s*linspace(-3,3);
p1=integral(f1s,-inf,inf,'Waypoints',w);
p2=integral(f1,-inf,inf,'Waypoints',w);
result=p2/p1;
  1 Comment
Paulo F.
Paulo F. on 20 Feb 2015
Thanks for your proposal, works beautifully.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

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

Start Hunting!