# randomdist.lua

## NAME

randomdist.lua - a few simple functions for generating random numbers.

## SYNOPSIS

 local R = require 'randomdist'
 grand1 = R.new_grand(10,3)
grand2 = R.new_grand(100,3)
for i = 1,20 do print( grand1(), grand2() ) end

gue_irand1 = R.new_gue_irand(4)
gue_irand2 = R.new_gue_irand(20)
for i = 1,20 do print( gue_irand1(), gue_irand2() ) end

for i = 1,20 do print(R.rayleigh_rand(3)) end

a = {'cold', 'cool', 'warm', 'hot'}
for i = 1,20 do print( R.randomget(a) ) end
for i = 1,20 do print(table.unpack( R.randomgetn(a, 2) )) end

a = {'a', 'b', 'c', 'd', 'e', 'f', 'g', 'h', }
zipf1 = R.new_zipf(a)
for i = 1,20 do print(R.zipf1()) end
zipf2 = R.new_zipf(8, 1.02)
for i = 1,20 do print(R.zipf2()) end

word2count = {
the=983, ['and']=421, of=340, to=286, I=263, it=252, -- etc
}
s, stddev = R.wordcount2zipf(word2count)
eo_words = {'la', 'kaj', 'de', 'al', 'mi', 'gxi'}
random_word = R.new_zipf(eo_words, s)
for i = 1,1000 do print(random_word()) end


## DESCRIPTION

This module implements in Lua a few simple functions for generating random numbers according to various distributions.

randomdist.lua is based on the PostScript module random.ps

## FUNCTIONS

new_grand (mean, stddev)

This function returns a closure, which is a function which you can then call to return a Gaussian (or Normal) Random distribution of numbers with the given mean and standard deviation.

It keeps some internal local state, but because it is a closure, you may run different Gaussian Random generators simultaneously, for example with different means and standard-deviations, without them interfering with each other.

It uses the algorithm given by Erik Carter which used to be at design.caltech.edu/erik/Misc/Gaussian.html

This algorithm generates results in pairs, but returns them one by one. Therefore if you are using math.randomseed to reset the random-number generator to a known state, and your code happens to make an odd number of calls to your closure, and you want your program to run consistently, then you should call your closure (eg: grand1) with the argument 'reset' each time you call math.randomseed. Eg:
   grand1 = R.new_grand(10,3)   ... grand1() ... etc ...   math.randomseed(244823040) ; grand1('reset') 

new_gue_irand (average)

This function returns a closure, which is a function which you can then call to return a Gaussian-Random-Ensemble distribution of integers.

The Gaussian Unitary Ensemble models Hamiltonians lacking time-reversal symmetry. Considering a hermitian matrix with gaussian-random values; from the ordered sequence of eigenvalues, one defines the normalized spacings
   s = (\lambda_{n+1}-\lambda_n) / <s>  where <s> = is the mean spacing. The probability distribution of spacings is approximately given by
   p_2(s) = (32 / pi^2) * s^2 * e^((-4/pi) * s^2)  These numerical constants are such that p_2 (s) is normalized: and the mean spacing is 1.
   \int_0^\infty ds p_2(s) = 1   \int_0^\infty ds s p_2(s) = 1 

Montgomery's pair correlation conjecture is a conjecture made by Hugh Montgomery (1973) that the pair correlation between pairs of zeros of the Riemann zeta function (normalized to have unit average spacing) is:
   1 - ({sin(pi u)}/{pi u}})^2 + \delta(u)  which, as Freeman Dyson pointed out to him, is the same as the pair correlation function of random Hermitian matrices.

rayleigh_rand (sigma)

This function returns a random number according to the Rayleigh Distribution, which is a continuous probability distribution for positive-valued random variables. It occurs, for example, when random complex numbers whose real and imaginary components are independent Gaussian distributions with equal variance and zero mean, in which case, the absolute value of the complex number is Rayleigh-distributed.
   f(x; sigma) = x exp(-x^2 / 2*sigma^2) / sigma^2 for x>=0  The algorithm contains no internal state, hence rayleigh_rand directly returns a number.

randomget (an_array)

This example gets a random element from the given array. For example, the following executes one of the four given functions at random:
   randomget( {bassclef, trebleclef, sharp, natural} ) () 

randomgetn (an_array, n)

This example returns an array containing n random elements, with distinct indices, from the given array.

new_zipf (an_array, s)
new_zipf (n, s)

This function returns a closure, which is a function which you can then call to return a Zipf-Distribution of array elements, or of integers.
The first example takes an array argument and returns a function which will return one of the items in the array, the first item being returned most frequently.
The second example takes an number argument and returns a function which will return a number from 1 to n, with 1 being the most frequent.
If s is not given it defaults to 1.0

s, stddev = wordcount2zipf (a_word_to_number_table)

This function can supply the s parameter used by new_zipf()
The argument is a table, for example:

city2population = {
Chongqing=30165500,
Shanghai=24183300,
Beijing=21707000,
Lagos=16060303,
Istanbul=15029231,
Karachi=14910352,    -- etc , etc ...
}

It returns two numbers: the Zipf-parameter s which best fits the data,
and the standard deviation stddev from which you can guess how reliable your parameter s is.

This module is available as a LuaRock in luarocks.org/modules/peterbillam so you should be able to install it with the command:

 \$ su
# luarocks install randomdist

or:

 # luarocks install http://www.pjb.com.au/comp/lua/randomdist-1.5-0.rockspec

The test script used during development is www.pjb.com.au/comp/lua/test_randomdist.lua

## CHANGES

 20180724 1.5 add wordcount2zipf
20171226 1.3 randomgetn allows n == array-size, meaning shuffle the array
20170707 1.1 grand('reset') more robust, and tested
20170706 1.0 first released version


## AUTHOR

Peter J Billam,   www.pjb.com.au/comp/contact.html

en.wikipedia.org/wiki/Normal_distribution
en.wikipedia.org/wiki/Random_matrix#Gaussian_ensembles
en.wikipedia.org/wiki/Random_matrix#Distribution_of_level_spacings
en.wikipedia.org/wiki/Montgomery%27s_pair_correlation_conjecture
{"ft_lang_label":"__label__en","ft_lang_prob":0.65150565,"math_prob":0.9779254,"size":6608,"snap":"2019-51-2020-05","text_gpt3_token_len":1805,"char_repetition_ratio":0.12113871,"word_repetition_ratio":0.075575024,"special_character_ratio":0.28208232,"punctuation_ratio":0.18329467,"nsfw_num_words":0,"has_unicode_error":false,"math_prob_llama3":0.9962602,"pos_list":[0,1,2],"im_url_duplicate_count":[null,null,null],"WARC_HEADER":"{\"WARC-Type\":\"response\",\"WARC-Date\":\"2020-01-23T06:05:50Z\",\"WARC-Record-ID\":\"<urn:uuid:1508afdb-cef2-4917-87cd-d83e142279b8>\",\"Content-Length\":\"12246\",\"Content-Type\":\"application/http; msgtype=response\",\"WARC-Warcinfo-ID\":\"<urn:uuid:2ed57a8d-7e6e-48fa-a88c-074911b5211f>\",\"WARC-Concurrent-To\":\"<urn:uuid:0130fa48-8c01-4798-be7d-8c7283caa22b>\",\"WARC-IP-Address\":\"173.255.248.22\",\"WARC-Target-URI\":\"https://pjb.com.au/comp/lua/randomdist.html\",\"WARC-Payload-Digest\":\"sha1:Y4UET46DWPGMBTI35Z3HCG5KF4EI7F3F\",\"WARC-Block-Digest\":\"sha1:AEWJFFETIP76UMJAXJGTUDHA7TNWORXA\",\"WARC-Identified-Payload-Type\":\"text/html\",\"warc_filename\":\"/cc_download/warc_2020/CC-MAIN-2020-05/CC-MAIN-2020-05_segments_1579250608295.52_warc_CC-MAIN-20200123041345-20200123070345-00353.warc.gz\"}"}