MathGroup Archive 2006

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: A conditional random number generation problem (please help me!)

  • To: mathgroup at smc.vnet.net
  • Subject: [mg65677] Re: A conditional random number generation problem (please help me!)
  • From: dh <dh at metrohm.ch>
  • Date: Fri, 14 Apr 2006 04:32:03 -0400 (EDT)
  • References: <e1ijuv$3oh$1@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

Hi,

if I understand correctly, you want a random number generator that gives 

numbers in [a,b] with a given density P.



This is a bit tricky. First note that to make life easier we can assume 

that the interval is [0,1]. This can then be mapped to [a,b] by:

a+x(b-a)



Further, Random[] gives random numbers in [0,1] uniformly distributed 

with density c (==constant). We are looking for a function F[x] that 

mappes the uniform distribution into one with the given density P[x]:

y=F[x]: [0,1] -> [0,1]

This must now be cast into mathematical notation. Assume we have a samll 

intervall dx that contains dn random numbers. This will be mapped into 

the interval:

dy== F'[x] dx

containing the same dn random numbers. Therefore,

dn== P[y] dy == P[F[x]] F'[x] dx == c dx

or

F'[x] == c/ P[ F[x] ]

That is the differential equation for F.

We also want that F mappes 0-> 0 and 1->1



Let's make an example:

we want a triangle distribution in [0,1]:



P[x_]:= 4 c If[x<0.5, x, 1-x]

the 4 c ensures that the numbers of random points stays constant. As 

this function has a kink at 0.5, we integrate from 0..0.5 and get the 

rest by symmetry. There is a further difficulty here, because P[0]==0, 

and 1/P[0] is indefinite. We circumvent this problem by replacing 0 by a 

very small number, say 10^-6.

Therefore, in [0,0.5] P[x]==x and c/P[F[x]]] == 0.25/F[x] and we can 

solve the differential equation:



fun0 = F /.  NDSolve[{F'[x] == 0.25/F[x], F[0] == 10^-6}, F,

                      {x, 0,0.5}][[1]];

fun[x_] := If[x < .5, fun0[x], 1 - fun0[1 - x]]

Plot[fun[x], {x, 0, 1}]



fun[x] will do the looked for mapping. We can test our result by simulation:



<< Graphics`Graphics`

d = Table[fun[Random[Real, {0, 1}]], {10000}];

Histogram[d]



Daniel







rjmachado3 wrote:

> I need to know the formula for the random function that return random

> numbers in a range of a and b integers [a,b] but that obey on a custom

> probability (possibly different!) for each integer number on this [a,b]

> range (of course the sum of all integer number probabilities are = 1!).

> Finally, what i want is the general function formula that simulate the

> random behavior (based on a custom probability value for each integer

> number between the [a,b] range. confuse? i hope not! please help me!!!!

> 

> what i know so far is that the function formula for generating a "pure"

> random number between [a,b] range is:

> 

> rand()*(b-a)+a

> 

> where rand() return a random number between 0 and 1.

> 



  • Prev by Date: 3D ContourGraphics alternative
  • Next by Date: Union, Sort at different levels of a list
  • Previous by thread: Re: A conditional random number generation problem (please help me!)
  • Next by thread: FindFit / NonlinearFit Problems