Introduction

Most nature processes are at random. From an infinite packaged order -the Big Bang- the Universe is relentlessly going to chaos, total absence of life, of movement. However randomness in lines, planes and volumes at a certain scale, seen from “upwards” may look like ordered, neat, with well defined lines surfaces and volumes. Dead tend to make all things the same. However, following a sequence of infinite “however”, we suppose that our Universe is causal and as such everything that happens in a system, either open or closed, is supposed to be “thermodynamically” deterministic, that is nature processes would be at last “pseudo random”, with all events being “at large” (at God scale perhaps?) predictable.

Pseudo Random Numbers: Randomness could be emulated by a computer program, for instance to generate (if binary) a sequence of (2^128 -1) strings of 128 bits each, all different but “at large” cycling in a predetermined sequence, from a given “seed” to go back to the same seed after (2^128-1) steps. They are predictable in the sense that given one string you may “guess” the string m steps ahead via a mathematic transformation, that is they are pseudo random generated by a virtual machines with (2^n – 1) different states S(j) at “time-steps” j’s. All those states are predictable from a given “seed” S(0) by an expression of the following type: S(j) = (T^j).S(0) where T^j is the “j power” of a unitary transformation T0 that obtain S(j+1) from S(j).

Random numbers are used for games, operations research, simulation, scientific experimenting and for the creation of crypto keys. In these cases, by obvious reasons, numbers should be unpredictable or at least extremely difficult to predict/unveil!. Ideally crypto keys should be created by True Random Numbers generators. True random numbers are those generated by nature process, for instance the position of a gel particle within a Brownian movement. Most nature processes are “entropy sources” in the sense that entropy is a state function.

Each state, for instance an ice cube, meanwhile its temperature is T has certain entropy that depends only of the temperature, no matter how the ice cube was formed. In the other extreme of complexity we may imagine a man of today with entropy that would depends on his “state” as a man alive, no matter how he achieves this state. Leaving this ice cube in an environment at a higher temperature, let’s say 20 Centigrade degrees, it will experiment a “degradation” becoming water, dying as ice.

If we isolate the ice cube and the room where the ice cube is on a table and consider both as a system we may measure that the air around the ice will frozen a little and we note that the ice cube start to melt becoming water that slides along the table wasting some energy: thermal energy transforms in kinetic energy. When a new equilibrium state is attained many things changed, no more ice cube being the most substantial change!. Globally we intuit that the order of the Universe, even though infinitesimally, has been degraded. This degradation is measured by the increase of its entropy, the more entropy the more the disorder is. In some extent entropy is also considered the “arrow of time” because as we go forward in time the entropy of an isolated system can only increase or remains the same, never decrease!. Entropy is in this concern the clock of the Universe.

Our concern is particularly related to Information Theory where entropy means how much randomness (or uncertainty) an event has. If we have to guess a number (or any object) out of a sequence of 16 our uncertainty “degree of ignorance” has a value of 16 if all numbers are equally probable. Now we are going to program a sort of “guessing game” with the aid of an informant that will provide us information -in the most economic way imagined- in order to reduce our ignorance gradually (we were tempted to say “bitwise”). First of all to impress the “public” we have to organize a little the scenario: the objects (letters) will be presented along a sequence in positions 0 to 15, representing pairs [position object] as depicted below:

Let’s program now the game supposing that we ignore the positional binary system. At most we were trained to know that a 0 means that something is in the lower half of a sequence and that a 1 means the contrary. To make things easier our game will deal with sets having 2^n objects, 2, 4, 8, 16, 32… 2^n, being n a natural number.

a) Someone of the public selects an object and informs our assistant either the object or the position chosen.

b) Our assistant must provide us the minimum but sufficient information to “guess” the chosen object, of course avoiding any unfair form.

Be the object chosen the number 11 (or the symbol &). Our assistant must provide us pieces of information to reduce our “ignorance” progressively; let’s say along a binary progression of the following type: 16 => 8 => 4 => 2 => 1. If now our assistant whisper the binary sequence [1010] in our ears to inform us that the right answer is the symbol located in position 10 in decimal we may say that we passed from an uncertainty measured by 16 to the certainty, namely from 16 to 1. The same information could be obtained progressively if previously we agreed with our assistant the procedure (as discussed above) : 1 meaning that the answer is in the upper half of the uncertainty range. So the sequence [1010] will mean:

Step 1: [1] the number is in the upper half [9 10 11 12 13 14 15 16] being 8 the level of our uncertainty after receiving the first “bit” of information.

Step2: [0] Once received the second bit as 0, according to the agreed procedure, our uncertainty is now reduced to 4, being certain that our number will be in the lower half, among numbers [9 10 11 12];

Step 3: [1] Our informant whisper now the third bit as a 1 and we know then that our number will be one of two, namely [11 12];

Step 4: [0] finally our informant whisper the fourth bit that permit us to be certain that the right answer is 11!.

Shannon defined Information Entropy as a variable related to the probabilities of events, outcomes of a random process like for instance flipping a coin. If all the symbol outcomes are equally likely then increasing the number of symbols should the entropy increases.

Information - Formal definition

Claude E. Shannon in its Mathematical Theory of Communications (a reprint from The Bell System Technical Journal) defines entropy H(x) in terms of a discrete random event x, with n possible states or outcomes as:

Where p(i) are probabilities of outcomes i’s. In our last examples n1=16, p1(i) = 1/16 for i=1, 2, 3,…, 16, if all outcomes were equally probable resulting H1=4, and if n2=2, (flipping a coin) p2(0-Head)= ½ , p2(1-Tail) = ½ , giving H1=1, being H in both cases the number of bits necessary to represent precisely the uncertainty of the guessing. Effectively with four bits we may represent 16 equally different outcomes, from 0000 to 1111 and with one bit the two flipping of coin outcomes, 0 for Head and 1 for Tail.

True Random Numbers

There are many Internet services providing True Random Numbers, using natural sources like radiation, Brownian movement, and quantum mechanics effects. Below we list some of them:

1. HotBits is a site from Switzerland at FermiLab. They said textually about themselves:

HotBits is an Internet resource that brings genuine random numbers, generated by a process fundamentally governed by the inherent uncertainty in the quantum mechanical laws of nature, directly to your computer in a variety of forms. HotBits are generated by timing successive pairs of radioactive decays detected by a Geiger-Müller tube interfaced to a computer.

The HotBits generation hardware produces data at a modest rate (about 30 bytes per second). Once the random bytes are delivered, they are immediately discarded--the same data will never be sent to any other user and no records are kept of the data at this or any other site.

A really good source of entropy is a radioactive source. The points in time at which a radioactive source decays are completely unpredictable, and can be sampled and fed into a computer, avoiding any buffering mechanisms in the operating system. In fact, this is what the HotBits people at Fermilab in Switzerland are doing.

Another sources of entropy could be atmospheric noise from a radio, like that used here at random.org , or even just background noise from an office or laboratory. The lavarand people at Silicon Graphics have been clever enough to use lava lamps to generate random numbers, so their entropy source not only gives them entropy, it also looks good! The latest random number generator to come online (both lavarand and HotBits precede random.org) is Damon Hart-Davis' Java EntropyPool which gathers random bits from a variety of sources including HotBits and random.org, but also from web page hits received by the EntropyPool's web server.

2. The LavaRND people at Silicon Graphics maintain a site that provides random numbers from lava lamps (kidding!). They proudly announce that LavaRND is a cryptographically sound random number generator. At its heart, they use not lava but some small thermal noise sources like some chips easy to find in Web cams like CCD.

3. Random.org . They use atmospheric noise to generate random numbers. They proudly announce themselves as the only service which offers a large (16K) block of numbers at once. They tuned a radio into a frequency where nobody is transmitting.

The atmospheric noise picked up by the receiver is fed into a workstation through the microphone port where it is sampled by a program at a given rate and the upper seven bits of each sample are discarded. The remaining bits are turned into a stream of bits with supposedly has a high entropy value. Skew Correction is performed on the bit stream, in order to ensure that there is an approximately even distribution of 0s and 1s.

4. ID Quantique . Photons - light particles - are sent one by one onto a semi-transparent mirror and detected. The exclusive events (reflection - transmission) are associated to "0" - "1" bit values. The operation of Quantis is continuously monitored to ensure immediate detection of a failure and disabling of the random bit stream. A High bit rate of 4Mbits/sec (up to 16Mbits/sec for PCI card) could be obtained.

The Quantum random bits generation process. See its justification

.

Random Numbers Products and Manufacturers

The Marsaglia CD-ROM

One of the best sources of random numbers is the George Marsaglia CD-ROM containing 600 megabytes of random numbers. These were produced by many sources (for example rap music).

Some suspected key ideas were used on it such as “if we have two collections of all possible bytes intentionally unordered (as octets) XOR them may kill all possible remaining order”. Marsaglia used three hardware sources identified by region origin: Canada, Germany and California.

Note: Some randomness gurus have found serious failures to this generator. One source of error could be in programs that manage data because some tricky details. For instance one of these experts found that Canadian and German generators failed when tested and inspecting the series they discover the following pattern: each byte containing a 10 was preceded by a byte containing a 13!. What could have happened?. One explanation was: when you write a program to write bytes to disk, you must open it as binary. If not the program thinks you are trying to write an end-of-line and the convention on a PC is to represent this by a carriage return (13) followed by a line feed (10).

Here are some manufacturers’ web sites:

• Colorado http://www.comscire.com/

• Holland http://valley.interact.nl/av/com/orion/home.html

• Canada http://www.tundra.com/ (find the data encryption section, then look under RBG1210. My device is an NM810 which is 2?8? RBG1210s on a PC card)

• Protego from Sweden, http://www.protego.se

• ID Quantique, http://www.idquantique.com/products/quantis.htm

Skew Correction

Recommended source: Request for Comments, RFC 1750

We are going to describe here a De-Skew technique, originally due to von Neumann. Suppose the probability of getting a 0 is equal to ½ plus a given error e and that all bits are independent. To eliminate skew, John von Neumann suggested to following algorithm:

• Read the bits two at a time.

• Skip 00’s and 11’s pairs.

• For 01 and 10 outputs the first bit.

Meaning to exam a bit stream as a sequence of non-overlapping pairs. We could then discard any 00 or 11 pairs found, interpret 01 as a 0 and 10 as a 1. Assume the probability of a 1 is 0.5 + e and the probability of a 0 is 0.5 - e where e is the eccentricity of the source. Then the probability of each pair is as follows:

So the probability of a 0 in the de-skewed sequence is the same as for a 1!. This technique will completely eliminate any bias but at the expense of taking an indeterminate number of input bits for any particular desired number of output bits. The probability of any particular pair being discarded is 0.5 + 2e^2 so the expected number of input bits to produce X output bits is X/(0.25 - e^2).

Defining Randomness

Sequences at random

Let’s consider infinite sequences of x’s pertaining to a GF2

.

Definition 1: Density D(x, n)

D(x, n) = 1/n. ∑ (x(i)), for i<n

Where the ∑ operator stands for Summa extended over i. This is like the “average” of ones. The limiting density of x’s is Lim(n)D(x, n) if that limit exists.

A true random sequence should have a limiting density ½ (though one might wonder why the limit should actually exist in the strict sense of convergence in mathematic analysis sense).

In fact, one would look for a certain degree of convergence: the D(x, n) should tend to ½

even though not too rapidly.

Definition 2: An infinite sequence of x’s pertaining to a GF2 is Mises random if the limiting density of any subsequence (xij) is ½ where the subsequence is selected by a rule named Auswahlregel.

For example, all the following sequences should have limiting density ½.

x0, x1, x2, . . . , xn, . . .

x0, x2, x4, . . . , x2n, . . .

x1, x4, x7, . . . , x3n+1, . . .

x0, x1, x4, . . . , xn2, . . .

Auswahlregeln Rule

This rule must be defined as a function N => N without any knowledge of x: otherwise

we can simply pick a subsequence of all 0’s. Now suppose we have a countable system of Auswahlregeln and our sequence passes all these tests. In other words, we have accounted for many increasing functions f : N => N and for each of these: lim(n)D((xf(n)), n) = ½ apply. Then we can think the sequence of x’s as being true random.

Some random numbers tables

1000 numbers ranged between 1 and 10^6 - First Half

1000 numbers ranged between 1 and 10^6 - Second Half

256 HX pairs

256 octal

Monte Carlo Techniques

Source: Department of Physics, Drexel University, Philadelphia .

This technique is fundamentally used to simulate statistic events and it derives from Integral Calculus. Let’s suppose we have to compute H by integrating function f over a given n-space as follows

H =.∫∫∫∫…f(x(1), x(2), x(3), …., x(n). dx(1).dx(2).dx(3)…..dx(n)

Where ∫’s stand for “Integrals” defined over domains in an n-dimensional space of the elementary computation (f.dv) where f stands for an n-dimensional function and dv for the infinitesimal hypercube defined by dv = dx(1).dx(2).dx(3)…..dx(n). In its discrete approximation It would mean m^n elementary computations of the core algorithm in a Cartesian space, namely m^n infinitesimal calculations of the form H <= H + f.dv along n loops of range m, where m stands for the number of sample points along each dimension supposedly all equal. If on the contrary for each dimension we have a particular m(i) the Cartesian space will have m(1).m(2).m(3)….m(n) points instead.

Note: That’s too much for m and n values easy to find in simulations, for example being n=6 and m=1000 for a 1000^6 = 10^18 grid for some weather forecast models!. We may navigate throughout this space systematically by “layers” at “brute force” mode (throughout all points) or following random walks. Monte Carlo approach estimates H via random samples within the Cartesian m^n space.

One dimensional case 1-D

For a given f defined over an interval [a b], the H computation will in pseudo code have the following form:

H=0, i=0

Do Core from i=0 till i=m, i++

Core: H = H + f(i)

H= delta.H

Where:

delta is the interval length divided in m parts, a practical infinitesimal.

f(i) is an array of m values along [a b]

Then the content of H at the end of computation will be

H = delta. ∑ (i) f(i).,

That reads delta multiplied by the ∑ of all array values. It’s the same as

H = delta.m.<f>,

Where <f> is the average of f along the interval [a b]. As delta.m is the length L of interval [a b] it results

H = L.<f>

So instead of computing H via “brute force” along m elementary computations, the Monte Carlo technique try to get a “good enough” statistical estimator of <f>. The economy rests on choosing a step larger than “delta” step. Let’s suppose that m = 10^8 points; applying the brute force method would imply to make 10^8 calculations an to keep someplace at hand equal number of f(i) values. What would be the resulting “error” in H if we use 1000 points instead of 10^8?: Perhaps not too much for our purpose. We may then take 1000 points spaced regularly along the interval [a b], equal number of corresponding f values and obtain <f> and its corresponding “variance” by the expression

Where <f^2> means the average of squares.

Pi Estimation by Monte Carlo

One of the classical Monte Carlo demos is the determination of number “pi” at a certain precision level. It’s based on the “quadrature” paradox of the antique world (Babylonians knew Pi with a reasonable precision as 25/8 = 3.125 that compared with Pi = 3.1416 gives us an error of 0.5%). If we integrate with a very elementary f (for instance f=1 within the X[-1, 1], Y[-1, 1] domain and zero outside it) along a circle of radio equal to 1, enclosed with a square of size equal to 2 we may find the following expression:

Area of the circle = (Pi)x1^2 = Pi

Area of the square = 4

Then area of the circle/area of the square = Pi/4, is a parameter that could be approached statistically by Monte Carlo techniques as follows. The idea is to “throw”/”pick” a zero dimension point chosen at random within the interval [-1, 1] in both dimensions X and Y of a Cartesian space. The probability p of these points to “fall” in either figure, the circle or the square is in relation to their areas

p(circle) = Pi/4

p(square) = 1

We may then generate N random pairs [x, y] of coordinates, for example between the range [-0.99999, 0.99999] and simply compute a Boolean variable P telling us if falling is inside or outside the circle. The pseudo code could be something like the one depicted below:

H=0

……………………….

x <= Random [-1 1]

y <= Random[-1 1]

If (x^2 + y^2>=1) P=1, and P=0 otherwise

H=H+P

…………………………

..

At the end of computation the estimation of Pi is given by four times H/N. As this is a pseudo Bernoulli process it will follow the Binomial Distribution that gives:

E(H) = E(<f>) = N.p

Variance = N.p.(1-p)

In numbers if we generate 10,000 pairs [x, y],

N=10,000,

E(H) = 10,000x3,1416/4 = 10,000x0.7854 = 7854

Totaling 7854 points within the circle with an expected standard deviation given by:

Variance = 10.000x(0.7854)x(0,2146) = 1685 => standard deviation = (1685)^ ½ = 41.

Statistics and Tests

NIST, http://csrc.nist.gov/rng/

The quality of random numbers can be measured in a variety of ways. One common method is to compute the Information Density, or entropy, in a series of numbers. The higher the entropy in a series of numbers is, the more difficult it is to predict a given number on the basis of the preceding numbers in the series. A sequence of good random numbers will have a high level of entropy, although a high level of entropy does not guarantee randomness. Paradoxically files compressed have a high level of entropy but as data is highly structured its randomness is not guaranteed being in fact considered not at random. Hence, for a thorough test of a random number generator, computing the level of entropy in the numbers alone is not enough.

Guide to Statistical Tests

http://csrc.nist.gov/rng/rng9.html

A total of sixteen statistical tests were developed, implemented and evaluated. The following describes each of the tests.

• Frequency (Monobits) Test

o The focus of the test is the proportion of zeroes and ones for the entire sequence. The purpose of this test is to determine whether that number of ones and zeros in a sequence are approximately the same as would be expected for a truly random sequence. The test assesses the closeness of the fraction of ones to ½, that is, the number of ones and zeroes in a sequence should be about the same.

• Test for Frequency within a Block

o The focus of the test is the proportion of zeroes and ones within M-bit blocks. The purpose of this test is to determine whether the frequency of ones is an M-bit block is approximately M/2.

• Runs Test

o The focus of this test is the total number of zero and one runs in the entire sequence, where a run is an uninterrupted sequence of identical bits. A run of length k means that a run consists of exactly k identical bits and is bounded before and after with a bit of the opposite value. The purpose of the runs test is to determine whether the number of runs of ones and zeros of various lengths is as expected for a random sequence. In particular, this test determines whether the oscillation between such substrings is too fast or too slow.

• Test for the Longest Run of Ones in a Block

o The focus of the test is the longest run of ones within M-bit blocks. The purpose of this test is to determine whether the length of the longest run of ones within the tested sequence is consistent with the length of the longest run of ones that would be expected in a random sequence. Note that an irregularity in the expected length of the longest run of ones implies that there is also an irregularity in the expected length of the longest run of zeroes. Long runs of zeroes were not evaluated separately due to a concern about statistical independence among the tests.

• Random Binary Matrix Rank Test

o The focus of the test is the rank of disjoint sub-matrices of the entire sequence. The purpose of this test is to check for linear dependence among fixed length substrings of the original sequence.

• Discrete Fourier Transform (Spectral) Test

o The focus of this test is the peak heights in the discrete Fast Fourier Transform. The purpose of this test is to detect periodic features (i.e., repetitive patterns that are near each other) in the tested sequence that would indicate a deviation from the assumption of randomness.

• Non-overlapping (Aperiodic) Template Matching Test

o The focus of this test is the number of occurrences of pre-defined target substrings. The purpose of this test is to reject sequences that exhibit too many occurrences of a given non-periodic (aperiodic) pattern. For this test and for the Overlapping Template Matching test, an m-bit window is used to search for a specific m-bit pattern. If the pattern is not found, the window slides one bit position. For this test, when the pattern is found, the window is reset to the bit after the found pattern, and the search resumes.

• Overlapping (Periodic) Template Matching Test

o The focus of this test is the number of pre-defined target substrings. The purpose of this test is to reject sequences that show deviations from the expected number of runs of ones of a given length. Note that when there is a deviation from the expected number of ones of a given length, there is also a deviation in the runs of zeroes. Runs of zeroes were not evaluated separately due to a concern about statistical independence among the tests. For this test and for the Non-overlapping Template Matching test, an m-bit window is used to search for a specific m-bit pattern. If the pattern is not found, the window slides one bit position. For this test, when the pattern is found, the window again slides one bit, and the search is resumed.

• Maurer's Universal Statistical Test

o The focus of this test is the number of bits between matching patterns. The purpose of the test is to detect whether or not the sequence can be significantly compressed without loss of information. An overly compressible sequence is considered to be non-random.

• Lempel-Ziv Complexity Test

o The focus of this test is the number of cumulatively distinct patterns (words) in the sequence. The purpose of the test is to determine how far the tested sequence can be compressed. The sequence is considered to be non-random if it can be significantly compressed. A random sequence will have a characteristic number of distinct patterns.

• Linear Complexity Test

o The focus of this test is the length of a generating feedback register. The purpose of this test is to determine whether or not the sequence is complex enough to be considered random. Random sequences are characterized by a longer feedback register. A short feedback register implies non-randomness.

• Serial Test

The focus of this test is the frequency of each and every overlapping m-bit pattern across the entire sequence. The purpose of this test is to determine whether the number of occurrences of the 2m m-bit overlapping patterns is approximately the same as would be expected for a random sequence. The pattern can overlap.

• Approximate Entropy Test

o The focus of this test is the frequency of each and every overlapping m-bit pattern. The purpose of the test is to compare the frequency of overlapping blocks of two consecutive/adjacent lengths (m and m+1) against the expected result for a random sequence.

• Cumulative Sum (Cusum) Test

o : The focus of this test is the maximal excursion (from zero) of the random walk defined by the cumulative sum of adjusted (-1, +1) digits in the sequence. The purpose of the test is to determine whether the cumulative sum of the partial sequences occurring in the tested sequence is too large or too small relative to the expected behavior of that cumulative sum for random sequences. This cumulative sum may be considered as a random walk. For a random sequence, the random walk should be near zero. For non-random sequences, the excursions of this random walk away from zero will be too large.

• Random Excursions Test

o : The focus of this test is the number of cycles having exactly K visits in a cumulative sum random walk. The cumulative sum random walk is found if partial sums of the (0,1) sequence are adjusted to (-1, +1). A random excursion of a random walk consists of a sequence of n steps of unit length taken at random that begin at and return to the origin. The purpose of this test is to determine if the number of visits to a state within a random walk exceeds what one would expect for a random sequence.

• Random Excursions Variant Test

o : The focus of this test is the number of times that a particular state occurs in a cumulative sum random walk. The purpose of this test is to detect deviations from the expected number of occurrences of various states in the random walk.