# Efficient binary vector

**URL:** <https://discourse.julialang.org/t/efficient-binary-vector/12066>\
**Category:** New to Julia\
**Tags:** question\
**Created:** [June 29, 2018, 2:06pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066 "2018-06-29T14:06:23Z")\
**Posts on this page:** 14\
**Page:** 2

<div class="post-metadata">

**Author:** ![foobar\_lv2](https://avatars.discourse-cdn.com/v4/letter/f/ee59a6/32.png) [@foobar\_lv2](https://discourse.julialang.org/u/foobar_lv2)\
**Post date:** [June 30, 2018, 5:02pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/21 "2018-06-30T17:02:12Z")

</div>

The random number generator is afaik already buffering its output. This is necessary because julia uses mersenne twister, which has an unholy large internal state.

So, as far as I gathered, you do _not_ need random access to you random vector? Then you should generate as you go, and possibly add your own buffer in front.

Secondly, how many random bits do you generate? Well, each symbol has an entropy of `-(p*log(p)+(p-1)*log(p-1))`, or in other words, 1e10 random coin-tosses with probability 1/3 : 2/3 only need 6.3e9 random bits. Generating 64 bits for each toin-coss is wasteful.

If you do need random access, I would also recommend to not store the random array. Instead store a seed, and use a PRNG with random access to the stream; for example `randomstream[n]=aes_encrpypt(key, n)`.

So, now your question should be: How do you transform a uniform random bitstream into a bernoulli-random bitstream? Preferably using only very few cycles / bit.

Unfortunately I don’t know. Something something asymmetric numeral systems should do the job; anyone else here knows?

* * *

AES as RNG:

My system takes (on julia 0.7) 0.84ns for mem-copying one word (8 byte), and 2.4ns for generating and storing one random word (via `Random.rand!` on large `Vector{UInt64}`) and ca 3ns for aes-encrypting one word (6ns ~ 12 cycles per encryption via `@btime`, without touching memory), using llvmcall and the aesni hardware instructions.

Amusingly this lines up almost perfectly with Agner Fog’s tables (11 aes instructions; these have a throughput of 1 and latency of 7 on Broadwell; this tells me that `@btime` is tight enough to fit into the reorder buffer, i.e. subsequent `@btime` iterations are actually done simultaneously. AES as a PRNG should probably do at least some buffering because we expect ~80 cycles latency for the encryption).

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [June 30, 2018, 5:41pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/22 "2018-06-30T17:41:07Z")

</div>

> [@foobar\_lv2](#):
>
> Generating 64 bits for each toin-coss is wasteful.

Possibly yes, but given that eg

```julia
sump(N, p) = sum(rand() < p for _ in 1:N)

@time sump(10^10, 1/3)

```

takes about 20s seconds on my (not very powerful) laptop, spending time on some sophisticated algorithm is also wasteful, in terms of developer time.

---

<div class="post-metadata">

**Author:** ![jlapeyre](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jlapeyre/32/4514_2.png) [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)\
**Post date:** [June 30, 2018, 6:39pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/24 "2018-06-30T18:39:17Z")

</div>

> Secondly, how many random bits do you generate? Well, each symbol has an entropy of -(p\*log§+(p-1)\*log(p-1)), or in other words, 1e10 random coin-tosses with probability 1/3 : 2/3 only need 6.3e9 random bits. Generating 64 bits for each toin-coss is wasteful.

Your expression for information entropy has a sign error and should have log2 not log (if it represents code, not math). Here is the correct formula.

```julia
julia> information_entropy(p) = -p*log2(p) - (1-p)*log2(1-p);
julia> information_entropy(1/2)
1.0
julia> information_entropy(1/3)
0.918295834054489

```

So a random string of Boolean samples with p=1/3 should be encodable with about 0.91 bits per bit in the string.

But, this is not the task. The task is to generate booleans with p=1/3 from a string with p=1/2.

(I deleted the previous post because Discourse sent it while I was trying to insert LaTeX)

---

<div class="post-metadata">

**Author:** ![semi-colon](https://avatars.discourse-cdn.com/v4/letter/s/43a26b/32.png) [@semi-colon](https://discourse.julialang.org/u/semi-colon)\
**Post date:** [June 30, 2018, 6:51pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/25 "2018-06-30T18:51:33Z")

</div>

Summing works for ordinary random walks, but not for reflected random walks. Also, I’m interested not in the final position of the walker, but rather in its maximum distance from the origin (over the entire duration of the walk). The code I have is:

```nohighlight
function rrw(N,p)
    X = BitArray(rand()<p for t in 1:N)    
    s = 0                                                                                
    m = 0  
    for t=1:N
      s = max(0,s+2*X[t]-1)        
      m = max(m,s)
    end
    return m
end

```

and suggestions for improving this would be welcome. The faster, the better. Thanks!

---

<div class="post-metadata">

**Author:** ![jlapeyre](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jlapeyre/32/4514_2.png) [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)\
**Post date:** [June 30, 2018, 7:03pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/26 "2018-06-30T19:03:01Z")

</div>

Using `X = [rand()<p for t in 1:N]` gives a `Vector{Bool}`. Your code runs about 10–20 % faster this way. Not a fantastic improvement.

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [June 30, 2018, 7:10pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/27 "2018-06-30T19:10:07Z")

</div>

```julia
function rrw2(N,p)
    s = 0                                                                                
    m = 0  
    for t=1:N
        X = rand() < p
        s = max(0,s+2*X-1)        
        m = max(m,s)
    end
    return m
end

```

then

```julia
julia> @time rrw(10^9, 1/3)
  4.117548 seconds (10 allocations: 119.214 MiB, 0.83% gc time)
27

julia> @time rrw2(10^9, 1/3)
  2.457705 seconds (6 allocations: 192 bytes)
28

```

---

<div class="post-metadata">

**Author:** ![benninkrs](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/benninkrs/32/3191_2.png) [@benninkrs](https://discourse.julialang.org/u/benninkrs)\
**Post date:** [June 30, 2018, 7:13pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/28 "2018-06-30T19:13:04Z")

</div>

This may not be relevant, but if p happens to be of the form k/2^n for integer k you could perform fast bitwise operations on n binary vectors to get the desired distribution. For example, for p = 1/4 take the AND of two uniformly distributed bit sequences.

---

<div class="post-metadata">

**Author:** ![jlapeyre](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jlapeyre/32/4514_2.png) [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)\
**Post date:** [June 30, 2018, 8:29pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/29 "2018-06-30T20:29:46Z")

</div>

This should probably be marked as the solution. The main question in the thread is “is it better or worse to generate an array of random numbers, rather than generating them as needed?” The answer is “worse”. At least in the vast majority of use cases.

---

<div class="post-metadata">

**Author:** ![foobar\_lv2](https://avatars.discourse-cdn.com/v4/letter/f/ee59a6/32.png) [@foobar\_lv2](https://discourse.julialang.org/u/foobar_lv2)\
**Post date:** [June 30, 2018, 9:06pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/30 "2018-06-30T21:06:50Z")

</div>

> [@jlapeyre](#):
>
> Your expression for information entropy has a sign error and should have log2 not log (if it represents code, not math).

Oops. Thanks.

> [@jlapeyre](#):
>
> But, this is not the task. The task is to generate booleans with p=1/3 from a string with p=1/2.

A p=1/2 random string decompresses into a p=1/3 random string 😉

Alas, Huffman needs to discretize probabilities to powers of 2, so one needs a fancy arithmetic/ANS/FSE decompressor. And I must admit that I can’t recall the details well enough to say whether this can be done chunkwise for the special case of two symbols. You could chunk 8 output bits, each symbol having inverse probability `3^k*1.5^(8-k)` , where k is the number of ones, and see whether some library permits setting the table explicitly and with high precision (e.g. each probability as `k/2^24`). Your generator would be something like `cat /dev/urandom | zstd --decompress --no-levzimpel --table=...` which should be pretty fast (the command names are made up; I know that zstd uses such a compression internally).

---

<div class="post-metadata">

**Author:** ![Juser](https://avatars.discourse-cdn.com/v4/letter/j/34f0e0/32.png) [@Juser](https://discourse.julialang.org/u/Juser)\
**Post date:** [July 1, 2018, 4:17am UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/31 "2018-07-01T04:17:56Z")

</div>

If you really care about speed, you can also probably parallelize your code by computing the max for different segments of the random walk and then connecting/aggregating the results at the end.

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [July 1, 2018, 4:35am UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/32 "2018-07-01T04:35:09Z")

</div>

How do you know where the walk will be in the different segments until you run the simulation? You can certainly parallelize over independent runs, though.

---

<div class="post-metadata">

**Author:** ![Juser](https://avatars.discourse-cdn.com/v4/letter/j/34f0e0/32.png) [@Juser](https://discourse.julialang.org/u/Juser)\
**Post date:** [July 1, 2018, 4:48am UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/33 "2018-07-01T04:48:07Z")

</div>

You can solve for the relative change from beginning to end of each segment as well as the max on each segment. You can then paste the endpoints together (the key step) to get the absolute levels of the walk at each of the endpoints and the absolute max on each segment. The absolute max over the whole walk is the absolute max over the segments.

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [July 1, 2018, 5:18am UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/34 "2018-07-01T05:18:46Z")

</div>

What I was trying to say is that you don’t know what the segments look like if eg. there’s a reflecting barrier; you only know if the system is translationally invariant and it doesn’t matter where you start.

---

<div class="post-metadata">

**Author:** ![Juser](https://avatars.discourse-cdn.com/v4/letter/j/34f0e0/32.png) [@Juser](https://discourse.julialang.org/u/Juser)\
**Post date:** [July 1, 2018, 4:50pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/35 "2018-07-01T16:50:38Z")

</div>

Good point. The presence of a reflecting barrier does make straightforward parallelism impossible.

[Previous page](https://discourse.julialang.org/t/efficient-binary-vector/12066.md?page=1)