# Efficiently resample a large vector?

**URL:** <https://discourse.julialang.org/t/efficiently-resample-a-large-vector/83177>\
**Category:** Performance\
**Created:** [June 22, 2022, 9:58am UTC](https://discourse.julialang.org/t/efficiently-resample-a-large-vector/83177 "2022-06-22T09:58:53Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![JedPhillips](https://avatars.discourse-cdn.com/v4/letter/j/aeb1de/32.png) [@JedPhillips](https://discourse.julialang.org/u/JedPhillips)\
**Post date:** [June 22, 2022, 9:58am UTC](https://discourse.julialang.org/t/efficiently-resample-a-large-vector/83177/1 "2022-06-22T09:58:53Z")

</div>

I’m doing some simulations and need to efficiently sample a set of items, and then subset/filter a large vector (1\>million entries) whose elements are contained in the sampled items. Here’s some minimal code.

```julia
Random.seed!(1234)
# generate some random data and corresponding element weights
DF = sample(1:15000, 1500000, replace = true)
wts = rand(Uniform(), 15000)
@time for s in 1:500
	# sample a subset of elements based on the weights
	samp = sample(1:15000, Weights(wts), 100, replace = false)
	# sample a single index from DF corresponding to each element
	sim_inds = [rand(findall(x-> x == i, DF)) for i in samp]
	# now use sim_inds for additional calculations....
end
 26.333749 seconds (788.80 k allocations: 8.999 GiB, 0.39% gc time, 0.39% compilation time)

```

I’m just suprised at how long this takes. I’ve tried various other approaches using dataframes, groupby, and so on, but the bottleneck always seems to be the findall step, or its equivalent, where I filter/match the sampled elements to the large vector. Any thoughts on how to speed this up?

EDIT: I’ve solved my own question in the act of posting it. StratifiedRandomSub in MLBase provides a massive speed boost:

```julia
@time begin
	# sample one index per element
	strat = hcat(collect(StratifiedRandomSub(DF, 15000, 500))...)
	# subsample the indices based on weights
	a = [sample(strat[:,i], Weights(wts), 100) for i in 1:500]
end
1.426467 seconds (7.62 M allocations: 903.555 MiB, 5.83% gc time, 2.54% compilation time)

```

But any further efficiency boost would be much appreciated!

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [June 22, 2022, 10:23am UTC](https://discourse.julialang.org/t/efficiently-resample-a-large-vector/83177/2 "2022-06-22T10:23:59Z")

</div>

You should never do benchmarking on non-`const` variables in global scope. Wrap your code in a function, and pass variables as input and output arguments.

Further performance tips here: [Performance Tips · The Julia Language](https://docs.julialang.org/en/v1/manual/performance-tips/) But the above tip is number 1.

And use the BenchmarkTools.jl package for benchmarking.

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [June 22, 2022, 10:29am UTC](https://discourse.julialang.org/t/efficiently-resample-a-large-vector/83177/3 "2022-06-22T10:29:13Z")

</div>

And one more thing: Try to make your example self-contained. For example, now `Random.seed!()`, `Uniform()` and `Weights()` are undefined and presumably come from some packages. Can you include the necessary code to load the correct libraries?

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [June 22, 2022, 2:45pm UTC](https://discourse.julialang.org/t/efficiently-resample-a-large-vector/83177/4 "2022-06-22T14:45:27Z")

</div>

A 1000x speedup by precomputing the indices :

```julia
using Random, Distributions, StatsBase, BenchmarkTools

# generate some random data and corresponding element weights
Random.seed!(1234)
DF = sample(1:15000, 1500000, replace = true)
wts = rand(Uniform(), 15000)
function original(DF,wts,N)
    for s in 1:N
        # sample a subset of elements based on the weights
        samp = sample(1:15000, Weights(wts), 100, replace = false)
        # sample a single index from DF corresponding to each element
        sim_inds = [rand(findall(x-> x == i, DF)) for i in samp]
        # now use sim_inds for additional calculations....
    end
end

w = Weights(wts)
indices = [findall(x -> x==i,DF) for i in 1:15000]
function my_try(indices,w,N)
    sim_inds = zeros(Int64,100)
    for s in 1:N
        samp = sample(1:15000, w, 100, replace = false)
        for i in 1:length(samp)
            sim_inds[i] = rand(indices[samp[i]])
        end
        # now use sim_inds for additional calculations....
    end
end 

@btime original(DF,wts,10) # 828ms
@btime my_try(indices,w,10) # 843μs

```

This will get even beter if you increase the number of repeats.

---

<div class="post-metadata">

**Author:** ![JedPhillips](https://avatars.discourse-cdn.com/v4/letter/j/aeb1de/32.png) [@JedPhillips](https://discourse.julialang.org/u/JedPhillips)\
**Post date:** [June 23, 2022, 10:40am UTC](https://discourse.julialang.org/t/efficiently-resample-a-large-vector/83177/5 "2022-06-23T10:40:53Z")

</div>

Thanks for the suggestion. But the timing ignores the time spent calculating the indices, which is still the bottleneck:

```julia
function get_inds(z::Vector{Int64}, r::UnitRange{Int64})
	indices = [findall(x -> x==i,z) for i in r]
	return indices
end

@btime get_inds(DF, 1:15000)
# 7.585 s (120040 allocations: 2.70 GiB)

```

---

<div class="post-metadata">

**Author:** ![rocco\_sprmnt21](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rocco_sprmnt21/32/20127_2.png) [@rocco\_sprmnt21](https://discourse.julialang.org/u/rocco_sprmnt21)\
**Post date:** [June 23, 2022, 10:52am UTC](https://discourse.julialang.org/t/efficiently-resample-a-large-vector/83177/6 "2022-06-23T10:52:19Z")

</div>

If this does what you are looking for, try a little …

```julia

d=Dict{Int, Vector{Int}}()
foreach(((i,e),)->push!(get!(()->Int[], d,e),i), enumerate(DF))
values(d) # or better simply d 

julia> Set(values(d))==Set(indices)
true

```