Non-Uniform Random Numbers in IronPython QuickStart Sample

Illustrates how to generate random numbers from a non-uniform distribution in IronPython.

This sample is also available in: C#, Visual Basic, F#.

Overview

This QuickStart sample demonstrates how to generate random numbers from non-uniform probability distributions using Numerics.NET.

The sample specifically illustrates the relationship between exponential and Poisson distributions in modeling events over time. It uses an exponential distribution to generate times between events, and then compares the resulting frequency of events per time unit with the theoretical Poisson distribution.

The code shows how to:

• Create and configure specific probability distributions (Exponential and Poisson)
• Use the MersenneTwister random number generator
• Generate random samples from a probability distribution
• Collect and analyze statistical data from the generated samples
• Compare empirical results with theoretical expectations

This practical example helps understand both the mechanics of generating non-uniform random numbers and how different probability distributions relate to each other in modeling real-world phenomena.

The code

import numerics

from System import Array

from Extreme.Statistics.Distributions import *
from Extreme.Statistics.Random import *

# Illustrates generating non-uniform random numbers
# using the classes in the Extreme.Statistics.Random
# namespace.

# Random number generators and the generation
# of uniform pseudo-random numbers are illustrated
# in the UniformRandomNumbers QuickStart Sample.

# In this sample, we will generate numbers from
# an exponential distribution, and compare summary
# results to what would be expected from 
# the corresponding Poisson distribution.

meanTimeBetweenEvents = 0.42

# We will use the exponential distribution to generate the time 
# between events. The number of events per unit time follows
# a Poisson distribution.

# The parameter of the exponential distribution is the time between events.
exponential = ExponentialDistribution(meanTimeBetweenEvents) 
# The parameter of the Poisson distribution is the mean number of events
# per unit time, which is the reciprocal of the time between events:
poisson = PoissonDistribution(1 / meanTimeBetweenEvents)

# We use a MersenneTwister to generate the random numbers:
random = MersenneTwister()

# The totals array will track the number of events per time unit.
totals = Array.CreateInstance(int, 15)

currentTime = 0
endOfCurrentTimeUnit = 1
eventsInUnit = 0

totalTime = 0
count = 0

while currentTime < 100000:
    timeBetween = exponential.Sample(random)
    totalTime += timeBetween 
    count = count + 1
    # Alternatively, we could have written
    #   timeBetween = random.NextDouble(exponential)
    # which would give an identical result.
    currentTime += timeBetween
    while currentTime > endOfCurrentTimeUnit:
        if eventsInUnit >= totals.Length:
            eventsInUnit = totals.Length-1
        totals[eventsInUnit] = totals[eventsInUnit] + 1
        eventsInUnit = 0
        endOfCurrentTimeUnit = endOfCurrentTimeUnit + 1
    eventsInUnit = eventsInUnit + 1

print "{0}", totalTime / count
# Now print the totals
print "# Events    Actual  Expected"
for i in range(0, totals.Length):
	expected = 100000 * poisson.Probability(i)
	print "{0:8}  {1:8}  {2:8.1f}".format(i, totals[i], expected)