# 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:** 20
**Page:** 1
### Author:  [@semi-colon](https://discourse.julialang.org/u/semi-colon)
#### Post date: [June 29, 2018, 2:06pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/1 "2018-06-29T14:06:23Z")
I wish to create a random vector X of length N with entries equal to either 1 or -1. For example,
```julia
using Distributions
X = 2.*rand(Bernoulli(p),N)-1
```
would work, where p=1/3. The catch, however, is that N=10^10. X, being Array{Int64}, consumes too much memory. But a trick like:
```julia
X = convert(Array{Int8}, 2.*rand(Bernoulli(p),N)-1)
```
also does not seem to help. Does anyone know a way to generate a Bernoulli(p) N-vector efficiently, i.e., as a binary array (not as a 64-byte integer array)? Thank you!
---
### Author:  [@dawbarton](https://discourse.julialang.org/u/dawbarton)
#### Post date: [June 29, 2018, 2:20pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/2 "2018-06-29T14:20:39Z")
How about `BitArray(rand() < p for x in 1:N)`?
---
### Author:  [@dawbarton](https://discourse.julialang.org/u/dawbarton)
#### Post date: [June 29, 2018, 2:23pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/3 "2018-06-29T14:23:43Z")
Or if you really do want the +1 or -1 entries, you can use `Int8[2*rand(Bernoulli(p))-1 for i in 1:N]` but that won’t be as space efficient.
---
### Author:  [@mbauman](https://discourse.julialang.org/u/mbauman)
#### Post date: [June 29, 2018, 2:45pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/4 "2018-06-29T14:45:53Z")
Definitely! `convert(Array{Int8}, 2.*rand(Bernoulli(p),N)-1)` probably isn’t helpful because it needs to first construct `2.*rand(Bernoulli(p),N)-1` before passing it to `convert`. Furthermore, `rand(Bernoulli(p),N)` is _also_ allocating a temporary before doing the multiplication. When you’re dealing with things this large, you want to be very careful about those temporaries. You can use broadcast fusion, but note that 10^10 bytes is still 10GB.
```julia
julia> using Compat
julia> function f(p, N)
A = Array{Int8}(undef, N)
A .= 2.*rand.(Bernoulli(p)) .- 1
return A
end
f (generic function with 1 method)
julia> @time f(.4, 10^7);
0.087543 seconds (7 allocations: 9.537 MiB)
julia> @time f(.4, 10^8);
0.840060 seconds (7 allocations: 95.368 MiB, 17.02% gc time)
julia> @time f(.4, 10^9);
7.077847 seconds (7 allocations: 953.675 MiB, 1.56% gc time)
```
It’s still gonna take some time for 10^10, and the vast majority of the time is being spent in `rand`. You can cut the space down by a factor of 8 if you use a `BitArray` like @dawbarton suggested, but the bigger savings come from using `bitrand` (which directly generates that `BitArray` with a p=0.5 distribution). A mapped array can make it still behave like an array of `Int` without any temporaries:
```julia
julia> using MappedArrays
julia> @time mappedarray(x->2*x-1, bitrand(10^10))
0.936435 seconds (1.07 k allocations: 1.164 GiB, 8.79% gc time)
10000000000-element MappedArrays.ReadonlyMappedArray{Int64,1,BitArray{1},##7#8}:
1
-1
…
```
---
### Author:  [@semi-colon](https://discourse.julialang.org/u/semi-colon)
#### Post date: [June 29, 2018, 3:22pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/5 "2018-06-29T15:22:47Z")
`bitrand` certainly is appealing, but I need to be able to set p=1/3 (or some value other than 1/2). Has anyone generalized `bitrand` in this manner? Grateful for your help…
---
### Author:  [@dawbarton](https://discourse.julialang.org/u/dawbarton)
#### Post date: [June 29, 2018, 4:15pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/6 "2018-06-29T16:15:47Z")
Looking at the implementation of `bitrand` it looks like it is tied to p=1/2. It generates a series of uniformally distributed `UInt64` values that are (essentially) reinterpreted as a binary string giving p=1/2. I personally can’t think of any way you can change that since most changes to the underlying random number generation will mean that the probability of each bit being true or false will change depending on the bit position.
My solution of `BitArray(rand() < p for x in 1:N)` is the same in terms of memory consumption but requires a call of the random number generator for each bit; this is the only (easy) way I can think of getting a variable probability but it is more costly computationally. In contrast, `bitrand` is exploiting the p=1/2 constraint to extract as much randomness from each call of the RNG as possible.
That said, if you only need to generate the BitVector once then that’s not so bad…
---
### Author:  [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)
#### Post date: [June 29, 2018, 4:16pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/7 "2018-06-29T16:16:13Z")
Consider an mmapped array, like this:
```julia
import Mmap
function random_bit_array(N, p)
_, io = mktemp()
A = Mmap.mmap(io, BitVector, N)
i = 0
@inbounds while i < N
ix = (i+1):min(i+8, N)
A[ix] .= rand(length(ix)) .< p
i += 8
end
A
end
@time A = random_bit_array(10^10, 1/3);
```
Will depend on your hard disk speed, for me this is around 100s.
---
### Author:  [@mbauman](https://discourse.julialang.org/u/mbauman)
#### Post date: [June 29, 2018, 5:01pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/8 "2018-06-29T17:01:46Z")
- `rand(Bernoulli(p),N)` is slow because it’s using 52 bits of randomness to generate each element.
- `bitrand(N)` is fast because it’s using 1 bit of randomness per element. Even better: it’s operating 64 bits at a time.
Of course, if you want `p!=0.5`, you’ll need more than one bit of randomness per element. But if you’re working with small fractions — particularly combinations of fractions with powers of two in the denominator — you can construct it yourself with a simple logical truth table:
```julia
bitrand(N) .& bitrand(N) # p = 1//4
bitrand(N) .& bitrand(N) .& bitrand(N) # p = 1//8
```
I imagine there’s some trick to efficiently get a 1//3, but my probability theory is too rusty to see it right now.
---
### Author:  [@Juser](https://discourse.julialang.org/u/Juser)
#### Post date: [June 29, 2018, 5:43pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/9 "2018-06-29T17:43:11Z")
Out of curiosity, in what context would one need to generate 10^10 bernoulli random variables? Also, if you really want to save time, you could generate, say 10^5 bernoulli random draws and then append that to itself 10^5 times. (Shuffle with each append if order really matters.) Technically, the distribution will not be a bernoulli, but it will be near enough to not matter in any application that I know of.
---
### Author:  [@semi-colon](https://discourse.julialang.org/u/semi-colon)
#### Post date: [June 30, 2018, 10:29am UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/10 "2018-06-30T10:29:33Z")
I’m interested in simulating 1D asymmetric simple random walks with reflection at the origin; see sections 1 & 2 of
> **[How Far Might We Walk at Random?](https://arxiv.org/abs/1802.04615)**
>
> This elementary treatment first summarizes extreme values of a Bernoulli random walk on the one-dimensional integer lattice over a finite discrete time interval. Both the symmetric (unbiased) and asymmetric (biased) cases are discussed. Asymptotic...
The low-memory approach (computing incremental steps one-by-one in sequence) is slow when N=10^10. If we instead efficiently store a vector X (computed before actually executing the walk), then a time savings should be possible.
---
### Author:  [@Per](https://discourse.julialang.org/u/Per)
#### Post date: [June 30, 2018, 10:52am UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/11 "2018-06-30T10:52:57Z")
Why would it be faster to first store a vector?
---
### Author:  [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)
#### Post date: [June 30, 2018, 11:19am UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/12 "2018-06-30T11:19:21Z")
> [@semi-colon](#):
>
> The low-memory approach (computing incremental steps one-by-one in sequence) is slow when N=10^10.
Are you sure about this? I see no intrinsic reason for this (unless, of course, the simulated quantity depends on the path globally and there is no online statistic). Perhaps some example code would make these things more concrete.
> [@semi-colon](#):
>
> instead efficiently store a vector X (computed before actually executing the walk)
Random number generation is deterministic once you set the seed, so it effectively “stores” the whole thing in a (relatively) small number of bytes.
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [June 30, 2018, 12:36pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/13 "2018-06-30T12:36:17Z")
> [@Per](#):
>
> Why would it be faster to first store a vector?
Maybe the poster is following habits from Python or Matlab, which train people into thinking that “vector” operations are always fast and loops are always slow. I agree that it seems unlikely to be the case here.
---
### Author:  [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)
#### Post date: [June 30, 2018, 1:49pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/14 "2018-06-30T13:49:33Z")
Maybe the RNG code is evicted from a cache if it is called only when a number is needed. But, I’ve never looked into this. Several years ago, I tested whether it was worth generating an array of random numbers for some Monte Carlo code for which the cpu time consumed by the RNG was significant. It made no difference. Still, I always used an array with the Mersenne Twister. I don’t recall if it was in the original code, or if I borrowed it from elsewhere.
EDIT: I definitely saw this method at least once in a prominent source. Maybe code from Knuth. You call a routine to get a sample. The routine uses a static array and an index and refills when necessary. The idea that this is the correct way may have propagated from there.
---
### Author:  [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)
#### Post date: [June 30, 2018, 2:02pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/15 "2018-06-30T14:02:16Z")
It’s easy to write routines that use 8/3 bits, on average, per sample (using rejection). I tried it last night and was not able to beat `rand() < 1/3` in speed. But, it still may be possible.
Maybe there is something analogous to the algorithm that generates two normally distributed samples.
---
### Author:  [@ScottPJones](https://discourse.julialang.org/u/ScottPJones)
#### Post date: [June 30, 2018, 2:02pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/16 "2018-06-30T14:02:38Z")
> [@jlapeyre](#):
>
> EDIT: I definitely saw this method at least once in a prominent source. Maybe code from Knuth. You call a routine to get a sample. The routine uses a static array and an index and refills when necessary. The idea that this is the correct way may have propagated from there.
That sort of technique is common when instruction & data caches come into play, and it can be faster in that case, greatly increasing locality of reference in both caches. If you have a large number of cores sharing L2/L3 caches on a physical CPU, and have a lot of tasks doing different things, it helps.
---
### Author:  [@DNF](https://discourse.julialang.org/u/DNF)
#### Post date: [June 30, 2018, 2:25pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/17 "2018-06-30T14:25:16Z")
> [@semi-colon](#):
>
> The low-memory approach (computing incremental steps one-by-one in sequence) is slow when N=10^10. If we instead efficiently store a vector X (computed before actually executing the walk), then a time savings should be possible.
Even if chunking the random generation were beneficial, surely you wouldn’t use such a gigantic cache that it would swallow that much memory? Why not use a reasonably sized cache that you periodically refill?
---
### Author:  [@Juser](https://discourse.julialang.org/u/Juser)
#### Post date: [June 30, 2018, 2:42pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/18 "2018-06-30T14:42:12Z")
@semi-colon, that sounds interesting. I agree with everyone else that there should be no performance improvement to pre-computing the N=10^10 random sequence. (Unless you’re pre-computing to ensure that you have repeated access to the identical sequence.)
---
### Author:  [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)
#### Post date: [June 30, 2018, 2:54pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/19 "2018-06-30T14:54:27Z")
Instructions are often stored in a separate cache. So its the RNG that should stay in cache.
EDIT: I just re-read. What I mean by “cache” is not a user allocated array, but cache in the CPU/memory system.
---
### Author:  [@ScottPJones](https://discourse.julialang.org/u/ScottPJones)
#### Post date: [June 30, 2018, 3:51pm UTC](https://discourse.julialang.org/t/efficient-binary-vector/12066/20 "2018-06-30T15:51:23Z")
If somebody is running a large number of workers all running the same code, using the same static tables, then it won’t matter so much, but if you have a mixed workload, then the RNG will tend to get pushed out of cache, which is why a (reasonably sized) buffer, can help performance (also when dealing with multiple threads).
We had the same issue for getting unique (sequential) ids, each process would get a batch (range) of ids, that it could hand out locally, instead of each process having to single thread every time a new id was required.
[Next page](https://discourse.julialang.org/t/efficient-binary-vector/12066.md?page=2)