algorithm - 将一组大整数转换为一组小整数

标签 algorithm combinatorics sampling random number-theory

我们如何重新编码一组严格递增(或严格递减)的正整数 P,以减少集合中整数之间可能出现的正整数的数量?

我们为什么要这样做:假设我们想对 P 进行随机抽样,但是 1.) P 太大而无法枚举,以及 2.) P 的成员以非随机方式相关,但以一种太复杂的方式无法采样。但是,当我们看到它时,我们知道它是 P 的成员。假设我们知道 P[0] 和 P[n],但不能接受枚举所有 P 或准确理解 P 的成员如何相关的想法。同样,出现在 P[0] 和 P[n] 之间的所有可能整数的数量比 P 的大小大很多倍,使得随机抽取 P 中成员的机会非常小。

示例:设 P[0] = 2101010101 & P[n] = 505050505。现在,也许我们只对 P[0] 和 P[n] 之间具有特定质量(例如,P[x] 中的所有整数总和为 Q 或更小,P 的每个成员的最大整数为 7 或更小)。所以,并非所有正整数 P[n] <= X <= P[0] 都属于 P。我感兴趣的 P 在下面的评论中讨论。

我尝试过的:如果 P 是一个严格递减的集合并且我们知道 P[0] 和 P[n],那么我们可以将每个成员视为从 P[ 0]。这样做会减少每个数字,也许会大大减少,并将每个成员保持为唯一的整数。对于我感兴趣的 P(下图),可以将 P 的每个减少值视为除以公分母 (9,11,99),这会减少 P 成员之间可能的整数数量。我已经发现结合使用,这些方法将所有 P[0] <= X <= P[n] 的集合减少几个数量级,从而有可能从所有正整数 P[n] 中随机抽取 P 的成员] <= X <= P[0] 仍然很小。

注意:应该清楚的是,我们必须了解一些关于 P 的信息。如果我们不了解,那基本上意味着我们不知道我们正在寻找什么。当我们随机抽取 P[0] 和 P[n] 之间的整数(重新编码或未重新编码)时,我们需要能够说“是的,它属于 P。”,如果确实如此。

一个好的答案可以大大增加我开发的计算算法的实际应用。评论 2 中给出了我感兴趣的 P 类型的示例。我坚持给予应有的信任。

最佳答案

虽然最初的问题是询问有关整数编码的非常通用的场景,但我认为不太可能存在完全通用的方法。例如,如果 P[i] 或多或少是随机的(从信息论的角度来看),如果有什么可行的话我会感到惊讶。

因此,让我们将注意力转向 OP 的实际问题,即生成恰好包含 K 个部分的整数 N 的分区。当将组合对象作为整数进行编码时,我们有必要尽可能多地保留组合结构。 为此,我们转向经典文本Combinatorial Algorithms by Nijenhuis and Wilf ,特别是第 13 章。事实上,在这一章中,他们展示了一个框架,用于从多个组合族中进行枚举和采样——包括 N 的分区,其中最大部分等于 K。使用分区之间众所周知的对偶性K部分和最大部分为K的分区(对Ferrers diagram进行转置),我们发现我们只需要对解码过程进行更改。

无论如何,这是一些源代码:

import sys
import random
import time

if len(sys.argv) < 4 :
    sys.stderr.write("Usage: {0} N K iter\n".format(sys.argv[0]))
    sys.stderr.write("\tN = number to be partitioned\n")
    sys.stderr.write("\tK = number of parts\n")
    sys.stderr.write("\titer = number of iterations (if iter=0, enumerate all partitions)\n")
    quit()

N = int(sys.argv[1])
K = int(sys.argv[2])
iters = int(sys.argv[3])

if (N < K) :
    sys.stderr.write("Error: N<K ({0}<{1})\n".format(N,K))
    quit()

# B[n][k] = number of partitions of n with largest part equal to k
B = [[0 for j in range(K+1)] for i in range(N+1)] 

def calc_B(n,k) :
    for j in xrange(1,k+1) :
        for m in xrange(j, n+1) :
            if j == 1 :
                B[m][j] = 1
            elif m - j > 0 :
                B[m][j] = B[m-1][j-1] + B[m-j][j]
            else :
                B[m][j] = B[m-1][j-1]

def generate(n,k,r=None) :
    path = []
    append = path.append

    # Invalid input
    if n < k or n == 0 or k == 0: 
        return []

    # Pick random number between 1 and B[n][k] if r is not specified
    if r == None :
        r = random.randrange(1,B[n][k]+1)

    # Construct path from r    
    while r > 0 :
        if n==1 and k== 1:
            append('N')
            r = 0   ### Finish loop
        elif r <= B[n-k][k] and B[n-k][k] > 0  : # East/West Move
            append('E')
            n = n-k
        else : #  Northeast/Southwest move
            append('N')
            r -= B[n-k][k]
            n = n-1
            k = k-1

    # Decode path into partition    
    partition = []
    l = 0
    d = 0    
    append = partition.append    
    for i in reversed(path) :
        if i == 'N' :
            if d > 0 : # apply East moves all at once
                for j in xrange(l) :
                    partition[j] += d
            d = 0  # reset East moves
            append(1) # apply North move
            l += 1            
        else :
            d += 1 # accumulate East moves    
    if d > 0 : # apply any remaining East moves
        for j in xrange(l) :
            partition[j] += d

    return partition


t = time.clock()
sys.stderr.write("Generating B table... ")    
calc_B(N, K)
sys.stderr.write("Done ({0} seconds)\n".format(time.clock()-t))

bmax = B[N][K]
Bits = 0
sys.stderr.write("B[{0}][{1}]: {2}\t".format(N,K,bmax))
while bmax > 1 :
    bmax //= 2
    Bits += 1
sys.stderr.write("Bits: {0}\n".format(Bits))

if iters == 0 : # enumerate all partitions
    for i in xrange(1,B[N][K]+1) :
        print i,"\t",generate(N,K,i)

else : # generate random partitions
    t=time.clock()
    for i in xrange(1,iters+1) :
        Q = generate(N,K)
        print Q
        if i%1000==0 :
            sys.stderr.write("{0} written ({1:.3f} seconds)\r".format(i,time.clock()-t))

    sys.stderr.write("{0} written ({1:.3f} seconds total) ({2:.3f} iterations per second)\n".format(i, time.clock()-t, float(i)/(time.clock()-t) if time.clock()-t else 0))

下面是一些性能示例(在 MacBook Pro 8.3、2GHz i7、4 GB、Mac OSX 10.6.3、Python 2.6.1 上):

mhum$ python part.py 20 5 10
Generating B table... Done (6.7e-05 seconds)
B[20][5]: 84    Bits: 6
[7, 6, 5, 1, 1]
[6, 6, 5, 2, 1]
[5, 5, 4, 3, 3]
[7, 4, 3, 3, 3]
[7, 5, 5, 2, 1]
[8, 6, 4, 1, 1]
[5, 4, 4, 4, 3]
[6, 5, 4, 3, 2]
[8, 6, 4, 1, 1]
[10, 4, 2, 2, 2]
10 written (0.000 seconds total) (37174.721 iterations per second)

mhum$ python part.py 20 5 1000000 > /dev/null
Generating B table... Done (5.9e-05 seconds)
B[20][5]: 84    Bits: 6
100000 written (2.013 seconds total) (49665.478 iterations per second)

mhum$ python part.py 200 25 100000 > /dev/null
Generating B table... Done (0.002296 seconds)
B[200][25]: 147151784574    Bits: 37
100000 written (8.342 seconds total) (11987.843 iterations per second)

mhum$ python part.py 3000 200 100000 > /dev/null
Generating B table... Done (0.313318 seconds)
B[3000][200]: 3297770929953648704695235165404132029244952980206369173   Bits: 181
100000 written (59.448 seconds total) (1682.135 iterations per second)

mhum$ python part.py 5000 2000 100000 > /dev/null
Generating B table... Done (4.829086 seconds)
B[5000][2000]: 496025142797537184410324290349759736884515893324969819660    Bits: 188
100000 written (255.328 seconds total) (391.653 iterations per second)

mhum$ python part-final2.py 20 3 0
Generating B table... Done (0.0 seconds)
B[20][3]: 33    Bits: 5
1   [7, 7, 6]
2   [8, 6, 6]
3   [8, 7, 5]
4   [9, 6, 5]
5   [10, 5, 5]
6   [8, 8, 4]
7   [9, 7, 4]
8   [10, 6, 4]
9   [11, 5, 4]
10  [12, 4, 4]
11  [9, 8, 3]
12  [10, 7, 3]
13  [11, 6, 3]
14  [12, 5, 3]
15  [13, 4, 3]
16  [14, 3, 3]
17  [9, 9, 2]
18  [10, 8, 2]
19  [11, 7, 2]
20  [12, 6, 2]
21  [13, 5, 2]
22  [14, 4, 2]
23  [15, 3, 2]
24  [16, 2, 2]
25  [10, 9, 1]
26  [11, 8, 1]
27  [12, 7, 1]
28  [13, 6, 1]
29  [14, 5, 1]
30  [15, 4, 1]
31  [16, 3, 1]
32  [17, 2, 1]
33  [18, 1, 1]

我会留给 OP 来验证此代码是否确实根据所需的(均匀)分布生成了分区。

编辑:添加了枚举功能的示例。

关于algorithm - 将一组大整数转换为一组小整数,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/14905190/

相关文章:

c - 使用模拟在一手 6 张牌中恰好有一张大牌

python-3.x - 如何在Python3中随机生成未观察到的数据

python - 带有事件持续时间的 Pandas TimeSeries

php - 使用 PHP 在 mySQL 中查找多个用户的帐户

在区域内查找三角形的算法

algorithm - 用于计算某个范围内的整数的数据结构?

R:组合学,一种组合的排列数

algorithm - 找到与数据一致的分数的所有可能组合

actionscript-3 - 麻烦让flashdevelop正确识别声音文件的采样率

algorithm - 计算if语句执行次数