下面是一个完成我所要求的脚本,就重新编码表示 N 和 K 部分的整数分区的整数而言。需要一种更好的重新编码方法才能使该方法适用于 K > 4。这绝对不是最佳或首选方法。然而,它在概念上很简单,很容易被认为是基本无偏见的。小 K 也很快。脚本在 Sage notebook 中运行良好,不调用 Sage 函数。它不是随机抽样的脚本。随机抽样本身不是问题。
方法:
1.) 将整数分区视为它们的和被连接在一起并根据第一个词法分区中最大和的大小用零填充,例如 [17,1,1,1] -> 17010101 & [5,5,5,5 ] -> 05050505
2.) 将得到的整数视为从最大整数中减去它们(即表示第一个词法分区的 int)。例如 17010101 - 5050505 = 11959596
3.) 将每个结果减少的整数视为除以公分母,例如 11959596/99 = 120804
所以,如果我们想选择一个随机分区,我们会:
1.) 选择一个介于 0 和 120,804 之间的数字(而不是介于 5,050,505 和 17,010,101 之间的数字)
2.) 将数字乘以 99 并从 17010101 中减去
3.) 根据我们如何将每个整数视为用 0 填充来拆分结果整数
优点和缺点:如问题正文中所述,这种特殊的重新编码方法不足以大大提高随机选择代表 P 成员的整数的机会。对于少量部分,例如 K < 5 并且基本上更大的总数,例如 N > 100,实现此概念的函数可能非常快,因为该方法避免了及时递归(蛇吃尾巴),这会减慢其他随机分区函数或使其他函数无法处理大 N。
在小 K 时,在考虑其余过程的速度时,绘制 P 成员的概率可能是合理的。再加上快速的随机抽取、解码和评估,此函数可以为 N&K 的组合(例如 N = 20000,K = 4)找到其他算法无法实现的统一随机分区。非常需要一种更好的方法来重新编码整数,以使其成为一种普遍强大的方法。
import random
import sys
首先,一些通常有用且简单明了的功能
def first_partition(N,K):
part = [N-K+1]
ones = [1]*(K-1)
part.extend(ones)
return part
def last_partition(N,K):
most_even = [int(floor(float(N)/float(K)))]*K
_remainder = int(N%K)
j = 0
while _remainder > 0:
most_even[j] += 1
_remainder -= 1
j += 1
return most_even
def first_part_nmax(N,K,Nmax):
part = [Nmax]
N -= Nmax
K -= 1
while N > 0:
Nmax = min(Nmax,N-K+1)
part.append(Nmax)
N -= Nmax
K -= 1
return part
#print first_partition(20,4)
#print last_partition(20,4)
#print first_part_nmax(20,4,12)
#sys.exit()
def portion(alist, indices):
return [alist[i:j] for i, j in zip([0]+indices, indices+[None])]
def next_restricted_part(part,N,K): # *find next partition matching N&K w/out recursion
if part == last_partition(N,K):return first_partition(N,K)
for i in enumerate(reversed(part)):
if i[1] - part[-1] > 1:
if i[0] == (K-1):
return first_part_nmax(N,K,(i[1]-1))
else:
parts = portion(part,[K-i[0]-1]) # split p
h1 = parts[0]
h2 = parts[1]
next = first_part_nmax(sum(h2),len(h2),(h2[0]-1))
return h1+next
""" *I don't know a math software that has this function and Nijenhuis and Wilf (1978)
don't give it (i.e. NEXPAR is not restricted by K). Apparently, folks often get the
next restricted part using recursion, which is unnecessary """
def int_to_list(i): # convert an int to a list w/out padding with 0'
return [int(x) for x in str(i)]
def int_to_list_fill(i,fill):# convert an int to a list and pad with 0's
return [x for x in str(i).zfill(fill)]
def list_to_int(l):# convert a list to an integer
return "".join(str(x) for x in l)
def part_to_int(part,fill):# convert an int to a partition of K parts
# and pad with the respective number of 0's
p_list = []
for p in part:
if len(int_to_list(p)) != fill:
l = int_to_list_fill(p,fill)
p = list_to_int(l)
p_list.append(p)
_int = list_to_int(p_list)
return _int
def int_to_part(num,fill,K): # convert an int to a partition of K parts
# and pad with the respective number of 0's
# This function isn't called by the script, but I thought I'd include
# it anyway because it would be used to recover the respective partition
_list = int_to_list(num)
if len(_list) != fill*K:
ct = fill*K - len(_list)
while ct > 0:
_list.insert(0,0)
ct -= 1
new_list1 = []
new_list2 = []
for i in _list:
new_list1.append(i)
if len(new_list1) == fill:
new_list2.append(new_list1)
new_list1 = []
part = []
for i in new_list2:
j = int(list_to_int(i))
part.append(j)
return part
最后,我们得到总 N 和部分数量 K。下面将按词汇顺序打印满足 N&K 的分区,以及相关的重新编码整数
N = 20
K = 4
print '#, partition, coded, _diff, smaller_diff'
first_part = first_partition(N,K) # first lexical partition for N&K
fill = len(int_to_list(max(first_part)))
# pad with zeros to 1.) ensure a strictly decreasing relationship w/in P,
# 2.) keep track of (encode/decode) partition summand values
first_num = part_to_int(first_part,fill)
last_part = last_partition(N,K)
last_num = part_to_int(last_part,fill)
print '1',first_part,first_num,'',0,' ',0
part = list(first_part)
ct = 1
while ct < 10:
part = next_restricted_part(part,N,K)
_num = part_to_int(part,fill)
_diff = int(first_num) - int(_num)
smaller_diff = (_diff/99)
ct+=1
print ct, part, _num,'',_diff,' ',smaller_diff
输出:
ct,分区,编码,_diff,smaller_diff
1 [17, 1, 1, 1] 17010101 0 0
2 [16, 2, 1, 1] 16020101 990000 10000
3 [15, 3, 1, 1] 15030101 1980000 20000
4 [15, 2, 2, 1] 15020201 1989900 20100
5 [14, 4, 1, 1] 14040101 2970000 30000
6 [14, 3, 2, 1] 14030201 2979900 30100
7 [14, 2, 2, 2] 14020202 2989899 30201
8 [13, 5, 1, 1] 13050101 3960000 40000
9 [13, 4, 2, 1] 13040201 3969900 40100
10 [13, 3, 3, 1] 13030301 3979800 40200
简而言之,最后一列中的整数可能要小得多。
为什么基于这种思想的随机抽样策略从根本上是无偏的:
N 的每个具有 K 个部分的整数分区对应于一个且仅一个重新编码的整数。也就是说,我们不会随机选择一个数字,对其进行解码,然后尝试重新排列元素以形成 N&K 的适当分区。因此,每个整数(无论是否对应于 N&K 的分区)都有相同的机会被绘制。目标是固有地减少不对应于 N 和 K 部分的分区的整数数量,从而使随机抽样的过程更快。