I wonder is there an algorithm, established package or some other way for inplace matrix sampling? Here is an example of what I am trying to achieve
using StableRNGs using StatsBase using BenchmarkTools rng = StableRNG(2020) a = rand(rng, 100, 500); a1 = copy(a) ids = sample(rng, axes(a, 2), size(a, 2));
It is easy to build new sampled matrix in the following manner
f1(a, ids) = a[:, ids]
but it allocates
julia> @btime f1($a, $ids); 39.015 μs (2 allocations: 390.70 KiB)
I can use
@views to get view of the sampled matrix, but unfortunately, later I need to change this matrix where each column should be treated independently, so
@view is not exactly what is needed.
On the other hand, simple assignment is too slow
function f2!(a, ids) b = @views a[:, ids] a .= b end julia> @btime f2!(x, $ids) setup=(x = copy(a)) evals = 1; 196.959 μs (3 allocations: 394.77 KiB)
and naive handwritten assignment is of course not working
function f3!(a, ids) b = @views a[:, ids] for j in axes(b, 2) for i in axes(b, 1) a[i, j] = b[i, j] end end end julia> b = f1(a, ids); julia> f3!(a, ids); julia> @assert a == b ERROR: AssertionError: a == b
It seems that it can be solved, but I can’t wrap my head around it.