posted April 28, 2005 02:46 PM            

quote:

---

Numbers generated by RND aren't really random, but are the result
of applying a pseudo-random transformation algorithm to a starting
("seed") value.

---

posted April 28, 2005 03:53 PM            

<<how the length or value of the Seed fed to the Randomizer impacts the period>>

The period is 2^32. The seed sets the point within that 2^32 numbers from which to start. The period remains 2^32 regardless of which point yo start.. unless you keep reseeding which is usually not a good idea.

Paul.

posted April 29, 2005 02:07 AM            

R1 = random number between 0 and 1 (for example)
R2 = random number between 0 and 0.5 (for example)
R3 = third random number, this one between 0 and 1 (again)
R_biased = R3*R1 + (1-R3)*R2

This skews the distribution of R1 towards the lower end of its range. You can influence the amount of skew by fiddling with the range of R2. If you set R2 to a constant, it will bias the distribution more significantly towards that number.

If you need more precise statistical control over the distribution - as opposed to my empirical "that looks OK for what I want" approach - I wouldn't recommend this solution!

------------------
Trevor

posted April 29, 2005 04:38 AM            

This is the (PB/DOS) source code to the Basically Speaking RNG article published in issue 07/1998.

Knuth

posted April 29, 2005 12:02 PM            

Getting various distributions of numbers isn't a job for the
random number generator itself, but rather you would want
an added function that transforms these numbers. For example,
here's a function that generates paired X,Y numbers that
are normally distributed, with any mean and standard deviation
you want, and any desired degree of correlation between X & Y.

The code is for PBCC but it should easily be translatable into
other versions of PB.

I don't have a similar function on hand for controlling skew,
but it should not be hard to find.

FUNCTION RandNormal (X() AS DOUBLE, Y() AS DOUBLE, PearsonR AS DOUBLE, MeanX AS DOUBLE, sdX AS DOUBLE, MeanY AS DOUBLE, sdY AS DOUBLE) AS LONG
      'Generate arrays of random normal numbers
      'from a population of Mean=0, sd=1 and any population Pearson r between -1 and +1
      'Note: You won't get **exactly** the r's, means, or sd's you ask for, but a sample
      'from a population
    LOCAL Pi2, BR, YY AS DOUBLE
    LOCAL LB,UB, I AS LONG
    LB=LBOUND (X): UB=UBOUND(X)
    IF PearsonR <-1 OR PearsonR >1 OR sdX<=0 OR sdY<=0 OR UB=LB OR UB<>UBOUND(Y) OR LB<>LBOUND(Y) THEN
        'illegal functions &/or math errors
        FUNCTION = 666    'error message
        EXIT FUNCTION
    END IF
    Pi2=ATN(1)*8

    BR = SQR(1 - PearsonR ^ 2)
    FOR I = LB TO UB
        X(I) = COS(Pi2 * RND) * SQR(-2 * LOG(RND))
         YY  = COS(Pi2 * RND) * SQR(-2 * LOG(RND))
        Y(I) = (YY * BR) + (PearsonR * x(I))
        'transform to requested population means * sd's
        X(I) = X(I)*sdX+MeanX  'a simple linear transform
        Y(I) = Y(I)*sdY+MeanY
    NEXT
    FUNCTION =0
END FUNCTION

posted April 29, 2005 02:34 PM            

I especially like the references in the postings and really
appreciate responses.

Emil,
Thanks for sharing your approach to biasing the distribution and for
reminding me of the Mersenne Twister.

Your RandNormal function isn't an approach I had thought about so
your idea should serve me well. No worries on it being PBCC. I do
all my Algorithm development in PBCC before moving it over to PBWin
so the approach is perfect for how I work.

Knuth,
Thanks for link. Is there any article text for that code snippet so
I can see if I might learn something?

Trevor,
I had tried using a distributed Bin approach that I sampled, but I
didn't like how it worked out. Your approach is a little more elegant
and may work when I try to do some Monte Carlo estimations.

Paul, thanks for the period on the generator.

posted April 29, 2005 05:45 PM            

This might be no more than what Trevor was talking about; I
wrote it before reading his comment. I don't, however, see it
as being either casual or makeshift, as one may have as many
data bins as one wishes, and as much precision in the
probability values as PowerBASIC can handle. But of course if
one wants any given precise degree of skewness one will have
to find a graph of a function that has that degree of skewness,
and work from there.

'RNDSKEW.BAS
'by Emil Menzel
'language: PBDOS
'actually this program could produce random numbers from
'any shape of theoretical distribution,
'but here the demo DATA statements
'are for a distribution with negative skew
'(for positive skew, reverse the order of the "numerical scores")
'For DATA enter first the number of data "bins"
'Then the numerical score for each bin... either high-to-low or vice versa
'then cumulative probabilities (CP), by bin... The final one should be 1.00
'and of course each CP must be  greater than (or equal to) the last one
$DIM ALL
DATA 9
DATA  10,  20,    30,   40,   50,   60,   70,   80,   90
DATA 0.02, 0.025, 0.05, 0.07, 0.09, 0.11, 0.15, 0.30, 1.00

DIM N as long, I as long, J as long
DIM Sum as double, expected as double, X as double
READ N
DIM CP(1 to N) as double
DIM Num(1 to N) as double
FOR I=1 TO N
    READ Num(I)
NEXT
FOR I=1 to N
    READ Cp(I)
NEXT
'==================================================================
'generate  a sample of N2 random numbers from the population with
'the distribution given in the DATA statements
DIM N2 as long
N2=1000
DIM Rn(1 to N2) as double
RANDOMIZE TIMER
FOR I=1 to N2
    X=RND  'a number between 0 and 1 (like a probability)
    'which bin does X fall into?
    FOR J=1 to N
        IF X<=CP(J) THEN  'if this bin
            Rn(I)=Num(J)  'transform the cumP value to the corresponding number
            EXIT FOR
        END IF
    NEXT J
NEXT I
'===============================================================
'check out these random data
DIM freq(N) as integer
'tally how many transformed RNDs fall in each bin
FOR I=1 to N
    X=Num(I)
    FOR J=1 to N2
       IF Rn(J)=X then INCR freq(I)
    NEXT
NEXT
'show the observed & expected frequencies
CLS
locate 12,1
print "score", "obs. freq.", "obs.cum.freq","'expected' cum.freq."
FOR I=1 to N
   Sum=Sum+freq(I)
   expected=CP(I)*N2
   PRINT Num(I),freq(I),Sum, expected
NEXT

------------------

posted April 29, 2005 08:30 PM            

This is similar to what I had tried. However, in my case
I used a sample distribution from a working process to create
a histogram of the output values. I then averaged the values
in the bins, some of which are negative, as the source values
that I used when the RNG value hit a bin window. In retrospect,
I probably should have used the median values, or picked from
the list of values in each bin with another random number.

In thinking about my approach, I now am beginning to think if I
had been more granular the results I needed might have looked
better than the jagged detail I ended up seeing. In simple
terms, my lack of finesse might be the reason why my bell's
tail was also so short.

Your example will be a good reflection to use when I compare it
to my approach.

Thank you for taking the time to make and send it.

posted April 29, 2005 11:56 PM            

I only wish it were easier reading & that the software was in
PB, not C.

posted April 30, 2005 10:23 AM            

Your link reference is a great resource. There are distributions on
that page that I had never seen before. I also found the entire
GSL Library and Table of Contents.

I hadn't seen this reference previously, so it should keep me busy
trying to figure out what is there.

Thanks for taking the time on this thread.

posted April 30, 2005 01:12 PM            

It was fun, not a chore.

If you learn how to link the libraries to PB, let me know.
It looks to me like they beat the classical "Numerical Recipes"
books for comprehensiveness.

No reply needed. I'm headed for the Gulf Shore very shortly.

Emil

posted May 02, 2005 05:07 AM            

quote:

---

Knuth,
Thanks for link. Is there any article text for that code snippet so
I can see if I might learn something?

---

The article was published on paper in Basically Speaking. I don't know if it is still possible to order back copies.

I translated it to German (with the permission of the author) and published it on my web site: http://www.softaware.de/rng.htm

There's a link to the original author (Jeff C. Jacobs) in that article, look for the link called "Email". Don't know if that address still works, though.

Two links that might be of interest, mentioned in the article:
http://random.mat.sbg.ac.at/literature/
http://www.iro.umontreal.ca/~lecuyer/

Knuth