Problem with extremely small real number

[Andrea]

Hi,
I need to calculate this probability P!/{n \choose p}, varying both n
and p, n takes values in this range [512:1024] and p in [2:12].
So i write this code in python:

def factorial(n):
result=1
if n==0: return 1
for i in xrange(1, abs(n)+1):
result = i*result
if n >= 0:
return result
else:
return -result
def binomial(n, k):
assert n>0 and isinstance(n, (int, long)) and isinstance(k,
(int,long))
if k < 0 or k > n:
return 0
if k == 0 or k == n:
return 1
return factorial(n) // (factorial(k) * factorial(n-k))

I want to call factorial(2)//binomial(1024,2) for example, in this way
trivially I obtain 0 as probability, how can I obtain the probability
forcing this division to output extremely small real numbers????

I want to produce a program that provide a set of data of this
probability on varying n and p, in order to plot a surface with this
data set with gnuplot, any comments or suggestion?
thanks in advance,
Andrea

 

[Duncan Smith]

Hi,
I need to calculate this probability P!/{n \choose p}, varying both n
and p, n takes values in this range [512:1024] and p in [2:12].
So i write this code in python:

def factorial(n):
result=1
if n==0: return 1
for i in xrange(1, abs(n)+1):
result = i*result
if n >= 0:
return result
else:
return -result
def binomial(n, k):
assert n>0 and isinstance(n, (int, long)) and isinstance(k,
(int,long))
if k < 0 or k > n:
return 0
if k == 0 or k == n:
return 1
return factorial(n) // (factorial(k) * factorial(n-k))

I want to call factorial(2)//binomial(1024,2) for example, in this way
trivially I obtain 0 as probability, how can I obtain the probability
forcing this division to output extremely small real numbers????

'//' is floor division.

Duncan

 

[Fabio Z Tessitore]

Il Mon, 03 Sep 2007 07:56:10 -0700, Andrea ha scritto:

[cut]

I want to call factorial(2)//binomial(1024,2) for example, in this way
trivially I obtain 0 as probability, how can I obtain the probability
forcing this division to output extremely small real numbers????

int to float will help

print float(factorial(2))/binomial(1024,2)

bye

 

[Arnaud Delobelle]

Hi,
I need to calculate this probability P!/{n \choose p}, varying both n
and p, n takes values in this range [512:1024] and p in [2:12].
So i write this code in python:

def factorial(n):
result=1
if n==0: return 1
for i in xrange(1, abs(n)+1):
result = i*result
if n >= 0:
return result
else:
return -result
def binomial(n, k):
assert n>0 and isinstance(n, (int, long)) and isinstance(k,
(int,long))
if k < 0 or k > n:
return 0
if k == 0 or k == n:
return 1
return factorial(n) // (factorial(k) * factorial(n-k))

I want to call factorial(2)//binomial(1024,2) for example, in this way
trivially I obtain 0 as probability, how can I obtain the probability
forcing this division to output extremely small real numbers????

I want to produce a program that provide a set of data of this
probability on varying n and p, in order to plot a surface with this
data set with gnuplot, any comments or suggestion?
thanks in advance,
Andrea

You're using integer division, that's why you probabilities are 0!
In python (up to 2.5 at list), if one integer/long is divided by
another, integer division is used.
But if one of the operands is a float, then the other one is converted
to a float and float division is used.

* convert your integers to floats first.
* Don't use //, use /
(// is 'true division', it will always give you the rounded answer)

[This can cause confusion, and I think it has been changed from v2.6
on - which I haven't got to try so correct me if I'm wrong]

Thus:
\3.8184261974584555e-06

BTW
p! / nCp = p! / (n! / p!(n-p)!) = (p!)^2 * (n! / (n-p)!) = (p!)^2 /
\prod_{i=n-p+1}^n i

Therefore:

def product(a, b):
return reduce(lambda x, y: x*y, xrange(a, b+1))

def f(n, p):
return float(product(1, p))**2 / product(n-p+1, n)
3.8184261974584555e-06

This will be more efficient.

 

[Roberto Bonvallet]

def binomial(n, k):
assert n>0 and isinstance(n, (int, long)) and isinstance(k,
(int,long))

Don't use assert to check whether the parameters have the right value.
assert should be used to claim that a certain condition always hold in
your program.
Prefer the following form, which will make debugging clearer:

if n <= 0:
raise ValueError("Wrong value of n")
if not isinstance(n, (int, long)):
raise TypeError("Wrong type of n")
if not isinstance(k, (int, long)):
raise TypeError("Wrong type of k")

Anyway, you shouldn't be checking the type of the parameters either.
It could be more useful if you cast to the desired type. This way,
your function will handle gracefully the case when you accidentally
pass a string as a parameter:

n = int(n)
k = int(k)
if n <= 0:
raise ValueError("Wrong value of n: %d" % n)

Best regards,