       • To: mathgroup at smc.vnet.net
• 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,

and 1/P 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 == 10^-6}, F,

{x, 0,0.5}][];

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

> 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

>

> 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