Skip to content

@pr4v33n /demo.py
Created

Embed URL

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP
cedric's challenge
"""
Challenge: http://beust.com/weblog/2012/02/16/a-new-coding-challenge-2/
Author: Praveen Kumar T
The Test:
--------
The probability distribution of the number of fixed points of a
uniformly distributed random permutation approaches a Poisson
distribution with expected value 1 as n grows.
he first n moments of this distribution are exactly those of
the Poisson distribution.
Source: http://en.wikipedia.org/wiki/Random_permutation#Statistics_on_random_permutations
"""
import scipy
import scipy.stats
def seq_producer(n):
"""
Produces a random shuffled array of [1..n]
"""
seq = scipy.arange(n)
scipy.random.shuffle(seq)
return seq
def count_fixed_points(n, seq):
"""
Count fixed points.
"""
orig_seq = scipy.arange(n)
return sum(p == q for p, q in zip(orig_seq, seq))
seq_size = 100 # sequence size to generate
sample_size = 100 # number of sequences to generate
# raw data of fixed points of the generated sequences
fp_data = [count_fixed_points(seq_size, seq_producer(seq_size)) for _ in range(sample_size)]
# number of moments to verify.
# Higher this value means more the sample data needed for accurate check
cmp_num_moments = 6
# central moments of the fixed point data
seq_moments = [scipy.stats.moment(fp_data, moment=i) for i in range(1, cmp_num_moments)]
# central moments of poisson distribution with mean = 1
# data courtesy: http://www.wolframalpha.com/input/?i=central+moments+of+poisson+distribution+with+mean+1
poisson_moments = [0, 1, 1, 4, 11, 41, 162, 715, 3425, 17722,
98253, 580317, 3633280, 24011157, 166888165,
1216070380, 9264071767, 73600798037,
608476008122, 5224266196935]
# sum the differences between corresponding moments
# Ideally this should be zero for a random sequence
print(sum(abs(x - y) for x, y in zip(seq_moments, poisson_moments)))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Something went wrong with that request. Please try again.