Julia 中更快的矢量比较

Faster vector comparison in Julia

我正在尝试使用 Julia 以最快的方式构造和比较两个长度相同的 01 随机向量,每个向量具有相同数量的零和一。

这都是为了对以下概率问题进行蒙特卡洛模拟

We have two independent urns, each one with n white balls and n black balls. Then we take a pair of balls, one of each urn, each time up to empty the urns. What is the probability that each pair have the same color?

我所做的是:

using Random
# Auxiliar function that compare the parity, element by element, of two 
# random vectors of length 2n
function comp(n::Int64)
  sum((shuffle!(Vector(1:2*n)) .+ shuffle!(Vector(1:2*n))).%2)
end

以上生成向量从 1 到 2n 的两个随机排列,逐个元素相加,对每个元素应用模 2,然后对剩余向量的所有值求和。然后我使用上面每个数字的奇偶性来模拟它的颜色:奇数黑色和白色偶数。

如果最终总和为零,则两个随机向量逐个元素具有相同的颜色。不同的结果表明这两个向量没有成对的颜色。

然后我设置了以下函数,它只是期望概率的蒙特卡洛模拟:

  # Here m is an optional argument that control the amount of random
  # experiments in the simulation
  function sim(n::Int64,m::Int64=24)
  # A counter for the valid cases
  x = 0
  for i in 1:2^m
    # A random pair of vectors is a valid case if they have the
    # the same parity element by element so
   if comp(n) == 0
    x += 1
   end
  end
  # The estimated value
   x/2^m
  end

现在我想知道是否有更快的方法来比较这些向量。我尝试了以下随机向量的替代构造和比较

shuffle!( repeat([0,1],n)) == shuffle!( repeat([0,1],n))

然后我相应地将代码更改为

comp(n)

经过这些更改,代码运行速度稍慢,这是我用函数 @time 测试的结果。我所做的其他更改是将 for 语句更改为 while 语句,但计算时间保持不变。

因为我不是程序员(实际上就在昨天我学习了一些 Julia 语言,并安装了 Juno 前端)那么可能是进行相同计算的更快方法。一些提示将不胜感激,因为蒙特卡洛模拟的有效性取决于随机实验的数量,因此计算速度越快,我们可以测试的值就越大。

为了更简洁,您可以直接生成 0/1 值的向量,然后让 Julia 检查向量是否相等,例如

function rndvec(n::Int64)
 shuffle!(vcat(zeros(Bool,n),ones(Bool,n)))
end

function sim0(n::Int64, m::Int64=24)
  sum(rndvec(n) == rndvec(n) for i in 1:2^m) / 2^m
end

正如 Bogumił Kamiński 所解释的那样,避免分配使代码更快(并且让 Julia 进行比较比他的代码更快)。

function sim1(n::Int64, m::Int64=24)
  vref = vcat(zeros(Bool,n),ones(Bool,n))
  vshuffled = vref[:]
  sum(shuffle!(vshuffled) == vref for i in 1:2^m) / 2^m
end

要更快地使用惰性计算和快速退出:如果第一个元素不同,您甚至不需要生成其余向量。 不过,这会使代码更加棘手。

我发现这有点不符合问题的精神,但你也可以做更多的数学运算。 生成了 binomial(2*n, n) 个可能的向量,因此您只需计算

function sim2(n::Int64, m::Int64=24)
   nvec = binomial(2*n, n)
   sum(rand(1:nvec) == 1 for i in 1:2^m) / 2^m
end

以下是我获得的一些时间:

@time show(("sim0", sim0(6, 21)))
@time show(("sim1", sim1(6, 21)))
@time show(("sim2", sim2(6, 21)))
@time test(("test", test(6, 2^21)))

("sim0", 0.0010724067687988281)  4.112159 seconds (12.68 M allocations: 1.131 GiB, 11.47% gc time)
("sim1", 0.0010781288146972656)  0.916075 seconds (19.87 k allocations: 1.092 MiB)
("sim2", 0.0010628700256347656)  0.249432 seconds (23.12 k allocations: 1.258 MiB)
("test", 0.0010166168212890625)  1.180781 seconds (2.14 M allocations: 98.634 MiB, 2.22% gc time)

这个问题的关键成本是 shuffle! 因此,为了最大限度地提高您可以使用的模拟速度(我将其添加为答案,因为它太长而无法评论):

function test(n,m)
  ref = [isodd(i) for i in 1:2n]
  sum(all(view(shuffle!(ref), 1:n)) for i in 1:m) / m
end

与其他答案中提出的代码有何不同:

  1. 您不必 shuffle! 两个向量; shuffle! 其中之一就足够了,因为比较的结果对于两个向量在独立洗牌后的任何相同排列都是不变的;因此我们可以假设一个向量在随机排列后重新排序以便它在第一个 n 条目中具有 trues 并且在最后一个 n 条目中具有 falses
  2. 我做 shuffle! in-place(即 ref 向量只分配一次)
  3. 我在向量的前半部分使用了all函数;这样,当我先点击 false 时,检查就会停止;如果我在第一个 n 条目中点击所有 true 我不必检查最后一个 n 条目因为我知道它们都是 false 所以我不必检查他们