Python 二项式系数

Question

import math
x = int(input("Enter a value for x: "))
y = int(input("Enter a value for y: "))

if y == 1 or y == x:
    print(1)

if y > x:
    print(0)        
else:
    a = math.factorial(x)
    b = math.factorial(y)
    div = a // (b*(x-y))
    print(div)  

这个二项式系数的程序是可以运行的，但当我输入两个相同的数字时，它应该等于1；而当y大于x时，它应该等于0。

Answer 1

使用递归定义是个不错的主意，就像Vadim Smolyakov的回答那样，再加上动态规划（DP）。对于动态规划，你可以使用functools模块里的lru_cache装饰器来帮助你：

import functools

@functools.lru_cache(maxsize = None)
def binom(n,k):
    if k == 0: return 1
    if n == k: return 1
    return binom(n-1,k-1)+binom(n-1,k)

Answer 2

所以，如果你搜索“在Python中实现二项式系数”，这个问题会首先出现。只有这个回答的第二部分包含了一个高效的实现方法，它依赖于乘法公式。这个公式能做到最少的乘法运算。下面的函数不依赖任何内置函数或导入库：

def fcomb0(n, k):
    '''
    Compute the number of ways to choose $k$ elements out of a pile of $n.$

    Use an iterative approach with the multiplicative formula:
    $$\frac{n!}{k!(n - k)!} =
    \frac{n(n - 1)\dots(n - k + 1)}{k(k-1)\dots(1)} =
    \prod_{i = 1}^{k}\frac{n + 1 - i}{i}$$

    Also rely on the symmetry: $C_n^k = C_n^{n - k},$ so the product can
    be calculated up to $\min(k, n - k).$

    :param n: the size of the pile of elements
    :param k: the number of elements to take from the pile
    :return: the number of ways to choose k elements out of a pile of n
    '''

    # When k out of sensible range, should probably throw an exception.
    # For compatibility with scipy.special.{comb, binom} returns 0 instead.
    if k < 0 or k > n:
        return 0

    if k == 0 or k == n:
        return 1

    total_ways = 1
    for i in range(min(k, n - k)):
        total_ways = total_ways * (n - i) // (i + 1)

    return total_ways

最后，如果你需要更大的数值，并且不介意牺牲一些准确性，斯特林近似可能是个不错的选择。

Answer 3

这是一个实际使用了正确公式的版本。 :)

#! /usr/bin/env python

''' Calculate binomial coefficient xCy = x! / (y! (x-y)!)
'''

from math import factorial as fac


def binomial(x, y):
    try:
        return fac(x) // fac(y) // fac(x - y)
    except ValueError:
        return 0


#Print Pascal's triangle to test binomial()
def pascal(m):
    for x in range(m + 1):
        print([binomial(x, y) for y in range(x + 1)])


def main():
    #input = raw_input
    x = int(input("Enter a value for x: "))
    y = int(input("Enter a value for y: "))
    print(binomial(x, y))


if __name__ == '__main__':
    #pascal(8)
    main()

...

这是我几年前写的一个不同版本的binomial()，它没有使用math.factorial()，因为在旧版本的Python中是没有这个函数的。不过，如果r不在0到n+1的范围内，它会返回1。

def binomial(n, r):
    ''' Binomial coefficient, nCr, aka the "choose" function 
        n! / (r! * (n - r)!)
    '''
    p = 1    
    for i in range(1, min(r, n - r) + 1):
        p *= n
        p //= i
        n -= 1
    return p

Answer 4

注意，从 Python 3.8 开始，标准库里新增了一个 math.comb 函数，用来计算二项式系数：

math.comb(n, k)

这个函数的意思是从 n 个物品中选择 k 个物品的方式有多少种，而且选择的时候不允许重复。它的计算公式是
n! / (k! (n - k)!)：

import math
math.comb(10, 5)  # 252
math.comb(10, 10) # 1

Answer 5

这个问题虽然比较老，但因为在搜索结果中排得很高，所以我想提一下，scipy有两个函数可以用来计算二项式系数：

1. scipy.special.binom()

2. scipy.special.comb()

   import scipy.special
   
   # the two give the same results 
   scipy.special.binom(10, 5)
   # 252.0
   scipy.special.comb(10, 5)
   # 252.0
   
   scipy.special.binom(300, 150)
   # 9.375970277281882e+88
   scipy.special.comb(300, 150)
   # 9.375970277281882e+88
   
   # ...but with `exact == True`
   scipy.special.comb(10, 5, exact=True)
   # 252
   scipy.special.comb(300, 150, exact=True)
   # 393759702772827452793193754439064084879232655700081358920472352712975170021839591675861424
   

需要注意的是，scipy.special.comb(exact=True)使用的是Python的整数，因此它可以处理非常大的结果！

在速度方面，这三个版本的表现会有些不同：

num = 300

%timeit [[scipy.special.binom(n, k) for k in range(n + 1)] for n in range(num)]
# 52.9 ms ± 107 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

%timeit [[scipy.special.comb(n, k) for k in range(n + 1)] for n in range(num)]
# 183 ms ± 814 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)each)

%timeit [[scipy.special.comb(n, k, exact=True) for k in range(n + 1)] for n in range(num)]
# 180 ms ± 649 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

（对于n = 300，二项式系数的值太大，无法用float64类型正确表示，如上所示）。