python - 有效地计算组合和排列

Question

我有一些代码来计算排列和组合，并且我正在努力让它更好地适用于大量数字。

我找到了一种更好的排列算法，可以避免大的中间结果，但我仍然认为我可以在组合方面做得更好。

到目前为止，我已经放入了一个特殊情况来反映 nCr 的对称性，但我仍然想找到一个更好的算法来避免调用 factorial(r)，这是一个不必要的大中间结果。如果没有这种优化，最后一个 doctest 会花费很长时间来计算阶乘（99000）。

任何人都可以提出一种更有效的计算组合的方法吗？

from math import factorial

def product(iterable):
    prod = 1
    for n in iterable:
        prod *= n
    return prod

def npr(n, r):
    """
    Calculate the number of ordered permutations of r items taken from a
    population of size n.

    >>> npr(3, 2)
    6
    >>> npr(100, 20)
    1303995018204712451095685346159820800000
    """
    assert 0 <= r <= n
    return product(range(n - r + 1, n + 1))

def ncr(n, r):
    """
    Calculate the number of unordered combinations of r items taken from a
    population of size n.

    >>> ncr(3, 2)
    3
    >>> ncr(100, 20)
    535983370403809682970
    >>> ncr(100000, 1000) == ncr(100000, 99000)
    True
    """
    assert 0 <= r <= n
    if r > n // 2:
        r = n - r
    return npr(n, r) // factorial(r)

score 27 · Accepted Answer

如果 n 离 r 不远，那么使用组合的递归定义可能会更好，因为 xC0 == 1 你只会有几次迭代：

这里相关的递归定义是：

nCr = (n-1)C(r-1) * n/r

这可以使用带有以下列表的尾递归很好地计算：

[(n - r, 0), (n - r + 1, 1), (n - r + 2, 2), ..., (n - 1, r - 1), (n, r)]

这当然很容易在 Python 中生成（我们省略了自 nC0 = 1 以来的第一个条目），izip(xrange(n - r + 1, n+1), xrange(1, r+1))请注意，这假设 r <= n 你需要检查它，如果不是，则交换它们。如果 r < n/2 则 r = n - r 也可以优化使用。

现在我们只需要使用带有reduce的尾递归来应用递归步骤。我们从 1 开始，因为 nC0 为 1，然后将当前值与列表中的下一个条目相乘，如下所示。

from itertools import izip

reduce(lambda x, y: x * y[0] / y[1], izip(xrange(n - r + 1, n+1), xrange(1, r+1)), 1)

score 20 · Accepted Answer

两个相当简单的建议：

为避免溢出，请在日志空间中执行所有操作。使用 log(a * b) = log(a) + log(b) 和 log(a / b) = log(a) - log(b) 的事实。这使得处理非常大的阶乘变得容易：log(n! / m!) = log(n!) - log(m!) 等。
使用 gamma 函数而不是阶乘。你可以在scipy.stats.loggamma. 这是一种比直接求和更有效的计算对数阶乘的方法。 loggamma(n) == log(factorial(n - 1))，同样，gamma(n) == factorial(n - 1)。

score 10 · Accepted Answer

在 scipy 中有一个尚未提及的功能：scipy.special.comb。根据您的 doctest 的一些快速计时结果（约 0.004 秒comb(100000, 1000, 1) == comb(100000, 99000, 1)），它似乎很有效。

[虽然这个具体问题似乎与算法有关，但问题是 python 中有一个数学 ncr 函数被标记为这个的副本......]

score 8 · Accepted Answer

8

如果您不需要纯 python 解决方案，gmpy2可能会有所帮助（gmpy2.comb非常快）。

于 2010-01-19T22:30:08.093 回答

score 6 · Accepted Answer

6

from scipy import misc
misc.comb(n, k)

应该允许您计算组合

于 2017-02-06T04:41:31.463 回答

score 5 · Accepted Answer

如果您正在计算 N 选择 K（我认为您正在使用 ncr 执行此操作），那么有一个动态编程解决方案可能会快很多。这将避免阶乘，而且如果您想以后使用，您可以保留该表。

这是它的教学链接：

http://www.csc.liv.ac.uk/~ped/teachadmin/algor/dyprog.html

不过，我不确定如何更好地解决您的第一个问题，抱歉。

编辑：这是模型。有一些非常搞笑的错误，所以它当然可以经受更多的清理。

import sys
n = int(sys.argv[1])+2#100
k = int(sys.argv[2])+1#20
table = [[0]*(n+2)]*(n+2)

for i in range(1,n):
    table[i][i] = 1
for i in range(1,n):
    for j in range(1,n-i):
        x = i+j
        if j == 1: table[x][j] = 1
        else: table[x][j] = table[x-1][j-1] + table[x-1][j]

print table[n][k]

score 5 · Accepted Answer

更有效的 nCr 解决方案 - 空间方面和精度方面。

中介 (res) 保证始终是 int 并且永远不会大于结果。空间复杂度为 O(1)（无列表、无拉链、无堆栈），时间复杂度为 O(r) - 恰好是 r 次乘法和 r 次除法。

def ncr(n, r):
    r = min(r, n-r)
    if r == 0: return 1
    res = 1
    for k in range(1,r+1):
        res = res*(n-k+1)/k
    return res

score 5 · Accepted Answer

对于 3.7 之前的 Python：

def prod(items, start=1):
    for item in items:
        start *= item
    return start


def perm(n, k):
    if not 0 <= k <= n:
        raise ValueError(
            'Values must be non-negative and n >= k in perm(n, k)')
    else:
        return prod(range(n - k + 1, n + 1))


def comb(n, k):
    if not 0 <= k <= n:
        raise ValueError(
            'Values must be non-negative and n >= k in comb(n, k)')
    else:
        k = k if k < n - k else n - k
        return prod(range(n - k + 1, n + 1)) // math.factorial(k)

对于 Python 3.8+：

有趣的是，组合功能的一些手动实现可能比math.comb()：

def math_comb(n, k):
    return math.comb(n, k)


def comb_perm(n, k):
    k = k if k < n - k else n - k
    return math.perm(n, k) // math.factorial(k)


def comb(n, k):
    k = k if k < n - k else n - k
    return prod(range(n - k + 1, n + 1)) // math.factorial(k)


def comb_other(n, k):
    k = k if k > n - k else n - k
    return prod(range(n - k + 1, n + 1)) // math.factorial(k)


def comb_reduce(n, k):
    k = k if k < n - k else n - k
    return functools.reduce(
        lambda x, y: x * y[0] // y[1],
        zip(range(n - k + 1, n + 1), range(1, k + 1)),
        1)


def comb_iter(n, k):
    k = k if k < n - k else n - k
    result = 1
    for i in range(1, k + 1):
        result = result * (n - i + 1) // i
    return result


def comb_iterdiv(n, k):
    k = k if k < n - k else n - k
    result = divider = 1
    for i in range(1, k + 1):
        result *= (n - i + 1)
        divider *= i
    return result // divider


def comb_fact(n, k):
    k = k if k < n - k else n - k
    return math.factorial(n) // math.factorial(n - k) // math.factorial(k)

所以实际上comb_perm()（用math.perm()and实现math.factorial()）实际上比math.comb()大多数时候都快。

请注意comb_reduce()，这很慢，与@wich的答案基本相同，而同样相对较慢，与@ZXX 的答案comb_iter()基本相同。

score 4 · Accepted Answer

如果您的问题不需要知道排列或组合的确切数量，那么您可以使用斯特林的阶乘近似值。

这将导致如下代码：

import math

def stirling(n):
    # http://en.wikipedia.org/wiki/Stirling%27s_approximation
    return math.sqrt(2*math.pi*n)*(n/math.e)**n

def npr(n,r):
    return (stirling(n)/stirling(n-r) if n>20 else
            math.factorial(n)/math.factorial(n-r))

def ncr(n,r):    
    return (stirling(n)/stirling(r)/stirling(n-r) if n>20 else
            math.factorial(n)/math.factorial(r)/math.factorial(n-r))

print(npr(3,2))
# 6
print(npr(100,20))
# 1.30426670868e+39
print(ncr(3,2))
# 3
print(ncr(100,20))
# 5.38333246453e+20

score 1 · Accepted Answer

from numpy import prod

def nCr(n,r):
    numerator = range(n, max(n-r,r),-1)
    denominator = range(1, min(n-r,r) +1,1)
    return int(prod(numerator)/prod(denominator))

score 0 · Accepted Answer

使用xrange()而不是range()会稍微加快速度，因为没有创建、填充、迭代和销毁中间列表。此外，reduce()与operator.mul.

score 0 · Accepted Answer

对于 N 选择 K，您可以使用帕斯卡三角形。基本上，您需要保留大小为 N 的数组来计算所有 N 选择 K 值。只需要添加。

score 0 · Accepted Answer

您可以输入两个整数并导入数学库以找到阶乘，然后应用 nCr 公式

import math
n,r=[int(_)for _ in raw_input().split()]
f=math.factorial
print f(n)/f(r)/f(n-r)

python - 有效地计算组合和排列

13 回答 13

Related

Reference