Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Using RAND in accordance with user defined probability distribution function?

Subject: Using RAND in accordance with user defined probability distribution function?

From: vinh

Date: 12 May, 2011 05:10:20

Message: 1 of 4

Hello,

I've done my fair amount of searching and the answer seems within reach, however my attempts at solving this problem haven't worked out too well.

I wish to use the rand function to generate random numbers w/ uniform distribution. My domain is from -1 to 1 so I know how to extend my result from rand to those boundaries. However, my distribution function is given as

g = 1 - x^2

My thought process has been to integrate my distribution function to get my density function x - x^3/3 ... the inverse of this should allow me to use my result from rand to generate random numbers in accordance with my distribution function g. However, running a histogram on the result is completely incorrect...I've also tried using the inverse of my distribution function g but that doesn't work either. Could someone help me with the correct way to think about this?

Subject: Using RAND in accordance with user defined probability distribution function?

From: pietro

Date: 12 May, 2011 06:52:05

Message: 2 of 4

"vinh" wrote in message <iqfq3s$gqa$1@newscl01ah.mathworks.com>...
> Hello,
>
> I've done my fair amount of searching and the answer seems within reach, however my attempts at solving this problem haven't worked out too well.
>
> I wish to use the rand function to generate random numbers w/ uniform distribution. My domain is from -1 to 1 so I know how to extend my result from rand to those boundaries. However, my distribution function is given as
>
> g = 1 - x^2
>
> My thought process has been to integrate my distribution function to get my density function x - x^3/3 ... the inverse of this should allow me to use my result from rand to generate random numbers in accordance with my distribution function g. However, running a histogram on the result is completely incorrect...I've also tried using the inverse of my distribution function g but that doesn't work either. Could someone help me with the correct way to think about this?

Hi,

I guess, this should help you
http://www.mathworks.se/matlabcentral/newsreader/view_thread/303811#822884

Bye

Pietro

Subject: Using RAND in accordance with user defined probability distribution function?

From: Roger Stafford

Date: 12 May, 2011 08:25:07

Message: 3 of 4

"vinh" wrote in message <iqfq3s$gqa$1@newscl01ah.mathworks.com>...
> Hello,
>
> I've done my fair amount of searching and the answer seems within reach, however my attempts at solving this problem haven't worked out too well.
>
> I wish to use the rand function to generate random numbers w/ uniform distribution. My domain is from -1 to 1 so I know how to extend my result from rand to those boundaries. However, my distribution function is given as
>
> g = 1 - x^2
>
> My thought process has been to integrate my distribution function to get my density function x - x^3/3 ... the inverse of this should allow me to use my result from rand to generate random numbers in accordance with my distribution function g. However, running a histogram on the result is completely incorrect...I've also tried using the inverse of my distribution function g but that doesn't work either. Could someone help me with the correct way to think about this?
- - - - - - - - - - - -
  I imagine you have already found that your g(x) needs to be multiplied by 3/4 to make it a valid probability density for x between -1 and +1, and you have undoubtedly found that the cumulative distribution function would then be

 G(x) = -1/4*x^3+3/4*x+1/2

  The hard part of your task is to then find the inverse to this function so that you can use 'rand'. In other words given y you are trying to find x such that y = -1/4*x^3+3/4*x+1/2 for any y in 0 <= y <= 1. Cubic equations like this do have explicit solutions but they are much messier than quadratic equations. The solution to this one for x in -1 <= x <= 1 and y in 0 <= y <= 1, courtesy of Cardan's formulas (1545), is as follows. Given a vector of y values, do this:

 a = 1-2*y;
 b = 2*sqrt(y.*(1-y));
 x = real((-1/2-1/2*sqrt(3)*i)*(a+b*i).^(1/3)+...
          (-1/2+1/2*sqrt(3)*i)*(a-b*i).^(1/3));

Then x will be the required corresponding roots. You can test this by computing -1/4*x.^3+3/4*x+1/2 to see if you come back to the original y values again (within round off error of course.)

  Hence if you do

 y = rand(1,n);

the above procedure will produce corresponding x values that possess the density distribution you desire, namely g(x) = 3/4*(1-x^2).

  Note: In theory the above values of x should already be real without applying the 'real' operator, but due to round off error they will usually have very small imaginary parts, so the 'real' operator is an easy way to remove these errors.

Roger Stafford

Subject: Using RAND in accordance with user defined probability distribution function?

From: Roger Stafford

Date: 12 May, 2011 15:44:20

Message: 4 of 4

"Roger Stafford" wrote in message <iqg5h3$el4$1@newscl01ah.mathworks.com>...
> .......
> a = 1-2*y;
> b = 2*sqrt(y.*(1-y));
> x = real((-1/2-1/2*sqrt(3)*i)*(a+b*i).^(1/3)+...
> (-1/2+1/2*sqrt(3)*i)*(a-b*i).^(1/3));
> .......
- - - - - - - - - -
  An alternative expression for that inverse function which avoids complex numbers using 'atan2' is:

 x = 2*cos((2*pi-atan2(2*sqrt(y.*(1-y)),1-2*y))/3);

  Another alternative, of course, is to use 'roots', always selecting the root that lies between -1 and +1.

Roger Stafford

Tags for this Thread

No tags are associated with this thread.

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us