高效的项目分类算法 (itertools/numpy)
我觉得这是一个常见的组合数学问题，但我找不到它的名字或相关资料。我现在在用Python和numpy来做这个，如果有快速的矩阵方法，我也能把它转过来。
简单来说，给定n个物品，我需要生成把它们放进m个箱子的所有方法。举个例子，把4个物品放进3个箱子，可能的组合会像这样[(4, 0, 0), (3, 1, 0), (3, 0, 1), (2, 2, 0), (2, 1, 1), ...]。这是一个总数固定的组合。
用itertools来实现这个是很简单的。
import itertools
def fixed_total_product(bins, num_items):
""" Return iterator of all item binning possibilities. """
return itertools.ifilter(lambda combo: sum(combo) == num_items,
itertools.product(xrange(num_items + 1), repeat=bins))
不过，我觉得后续在循环中进行计算会很低效。用2D的numpy数组来处理会快很多，但我不知道怎么高效地构建这个数组。我可以遍历ifilter的结果，生成可能性的列表，然后用这个来构建数组，但这看起来很浪费。
我猜最好的办法是用“numpy的方式”来构建一切，但我不太确定怎么做。在stackoverflow上有一个快速的乘积实现：使用numpy构建两个数组的所有组合的数组。我猜可以修改这个方法，只输出和正确的乘积。数组的大小应该是((m-1) + n)选择n，因为有m-1个箱子的边界。
有什么想法吗？如果有基准测试那就更好了，但不是必须的。
3 个回答
我知道我在这里提到的是一个比较老的话题，但我希望能有一个好的解决方案，哪怕现在发帖的人不再关注。
这个问题其实就像是在代码中表示多变量多项式。比如说，想象一下有4个物品要放进3个箱子里，你可以把这个问题转化为展开 (x+y+z)^4 得到的多项式，结果大概是 x^4*y^0*z^0 + 4*x^3*y^1*z^0 + ...，这里的x、y、z的指数就是数组中的元素：
[4,0,0],
[3,1,0],
[3,0,1],
[2,2,0],
[2,1,1],
[2,0,2],
[1,3,0],
[
如果你仔细观察，就会发现，比如说在 3 右边的元素是 (y+z)^1 的指数，在 2 右边的元素是 (y+z)^2 的指数，依此类推。更一般地说，在 p 右边的就是 (y+z)^(4-p) 的指数。这表明有一种递归的结构：假设 e(m,p) 表示 m 个变量（也就是箱子）和 p 个总指数（也就是物品），那么你可以通过选择第一个指数 q 从 0 到 p，剩下的 m-1 个变量的指数就可以用 e(m-1,p-q) 来表示。
在 python 中，使用 numpy 可以这样来实现：
def multiindex_exact(m, p):
if m == 0:
return np.zeros((1 if p==0 else 0, 0), np.int8)
else:
I = np.zeros((0, m), np.int8)
for q in reversed(range(0, p + 1)):
J = multiindex_exact(m - 1, p - q)
Jn = np.full((J.shape[0], 1), q)
I = np.vstack((I, np.hstack((Jn, J))))
return I
当然，你可以通过预先计算数组的大小来提高效率，并直接填入值。
如果你需要分配的物品不是严格的，而是 最多 有 p 个物品放入 m 个箱子中，使得所有箱子的总和始终 <=p，你可以使用以下类似的代码。
def multiindex_lower(m, p):
if m == 0:
return np.zeros((1, 0), np.int8)
else:
I = np.zeros((0, m), np.int8)
for q in range(0, p + 1):
J = multiindex_lower(m - 1, q)
Jn = q - J.sum(1).reshape((J.shape[0], 1))
I = np.vstack((I, np.hstack((J, Jn))))
return I
其实用numpy有一些更快的方法，可以用几个小技巧来实现。你可以从numpy.indices开始。这个函数基本上和itertools.product是一样的，只不过你需要把它和rollaxis结合起来使用。Sven Marnach在这个问题中的回答是个很好的例子（不过他最后一个例子有个小错误，你要用的应该是numpy.indices((len(some_list) + 1), * some_length...）
不过，对于这种情况，使用itertools可能会更容易理解。
numpy.fromiter比直接把东西转换成列表要快一点，特别是当你告诉它迭代器中有多少个项目时。它的主要优点是占用的内存少，这在处理大迭代器时非常有用。
下面是一些例子，展示了如何使用numpy.indices这个技巧，以及不同的方法把迭代器转换成numpy数组：
import itertools
import numpy as np
import scipy.special
def fixed_total_product(bins, num_items):
return itertools.ifilter(lambda combo: sum(combo) == num_items,
itertools.product(xrange(num_items + 1), repeat=bins))
def fixed_total_product_fromiter(bins, num_items):
size = scipy.special.binom(bins - 1 + num_items, num_items)
combinations = fixed_total_product(bins, num_items)
indicies = (idx for row in combinations for idx in row)
arr = np.fromiter(indicies, count=bins * int(size), dtype=np.int)
return arr.reshape(-1, bins)
def fixed_total_product_fromlist(bins, num_items):
return np.array(list(fixed_total_product(bins, num_items)), dtype=np.int)
def fixed_total_product_numpy(bins, num_items):
arr = np.rollaxis(np.indices((num_items + 1,) * bins), 0, bins + 1)
arr = arr.reshape(-1, bins)
arr = np.arange(num_items + 1)[arr]
sums = arr.sum(axis=1)
return arr[sums == num_items]
m, n = 5, 20
if __name__ == '__main__':
import timeit
list_time = timeit.timeit('fixed_total_product_fromlist(m, n)',
setup='from __main__ import fixed_total_product_fromlist, m, n',
number=1)
iter_time = timeit.timeit('fixed_total_product_fromiter(m, n)',
setup='from __main__ import fixed_total_product_fromiter, m, n',
number=1)
numpy_time = timeit.timeit('fixed_total_product_numpy(m, n)',
setup='from __main__ import fixed_total_product_numpy, m, n',
number=1)
print 'All combinations of {0} items drawn from a set of {1} items...'.format(m,n)
print ' Converting to a list and then an array: {0} sec'.format(list_time)
print ' Using fromiter: {0} sec'.format(iter_time)
print ' Using numpy.indices: {0} sec'.format(numpy_time)
关于时间方面：
All combinations of 5 items drawn from a set of 20 items...
Converting to a list and then an array: 2.75901389122 sec
Using fromiter: 2.10619592667 sec
Using numpy.indices: 1.44955015182 sec
你会发现这些方法都不是特别快。
你使用的是一种暴力算法（生成所有可能的组合，然后再筛选），我在这里的numpy示例也是这样做的。
其实可能有更高效的方法来实现这个！不过我并不是计算机科学或数学方面的人，所以我不知道有没有一种众所周知的算法可以在不生成所有可能组合的情况下做到这一点……
无论如何，祝你好运！
这个内容是关于一个数学公式的，它用递归的方法来计算组合数。组合数C(n, k)表示从n个物品中选出k个物品的不同方式。公式的意思是：如果你想知道从n个物品中选k个的方式，可以先计算从n-1个物品中选k个的方式，再加上从n-1个物品中选k-1个的方式。
这里还提到了一种叫做“备忘录”的技术，这种技术可以帮助我们记住已经计算过的结果，以便下次用的时候不需要再计算，节省时间。
另外，使用了numpy这个库，它是Python中一个非常流行的数学计算库，能让我们更高效地进行数值运算。
最后，文中提到的“时间（秒）”是指计算组合数C(20, 5)所花费的时间。
import numpy as np
def binnings(n, k, cache={}):
if n == 0:
return np.zeros((1, k))
if k == 0:
return np.empty((0, 0))
args = (n, k)
if args in cache:
return cache[args]
a = binnings(n - 1, k, cache)
a1 = a + (np.arange(k) == 0)
b = binnings(n, k - 1, cache)
b1 = np.hstack((np.zeros((b.shape[0], 1)), b))
result = np.vstack((a1, b1))
cache[args] = result
return result
if __name__ == '__main__':
import timeit
print timeit.timeit('binnings(20, 5, {})', setup='from __main__ import binnings', number=1)
对于(20, 5)的计算，所需的时间（秒）：
0.0129251480103