# Quantile! much faster when called twice with a scalar than once with a vector!

**URL:** <https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188>\
**Category:** General Usage\
**Tags:** statistics\
**Created:** [July 23, 2021, 3:03pm UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188 "2021-07-23T15:03:09Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [July 23, 2021, 3:03pm UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188/1 "2021-07-23T15:03:09Z")

</div>

```julia
x = rand(Float32, 10000000);

@time quantile!(x, [0.02, 0.98])
  0.196571 seconds (2 allocations: 160 bytes)
2-element Vector{Float64}:
 0.02001022338867188
 0.9800319683551788

@time quantile!(x, 0.02)
  0.022732 seconds (1 allocation: 16 bytes)
0.02001022338867188

@time quantile!(x, 0.98)
  0.022851 seconds (1 allocation: 16 bytes)
0.9800319683551788

```

---

<div class="post-metadata">

**Author:** ![stillyslalom](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stillyslalom/32/45687_2.png) [@stillyslalom](https://discourse.julialang.org/u/stillyslalom)\
**Post date:** [July 23, 2021, 3:11pm UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188/2 "2021-07-23T15:11:27Z")

</div>

`quantile!` modifies the input, so your second two cases don’t need to spend any time sorting the input vector.

```julia
julia> v = [3, 2, 1]
3-element Vector{Int64}:
 3
 2
 1

julia> quantile!(v, 0.5)
2.0

julia> v
3-element Vector{Int64}:
 1
 2
 3

```

That said, I can confirm that performance is much worse than it ought to be:

```julia
julia> @btime quantile($x, [0.02, 0.98])
  2.145 s (4 allocations: 38.15 MiB)

julia> @btime quantile($x, 0.02)
  99.073 ms (2 allocations: 38.15 MiB)

```

---

<div class="post-metadata">

**Author:** ![Salmon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/salmon/32/22968_2.png) [@Salmon](https://discourse.julialang.org/u/Salmon)\
**Post date:** [July 23, 2021, 3:53pm UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188/3 "2021-07-23T15:53:30Z")

</div>

To add to this: I am not super familiar with quantile!, but this is not due to the vector possibly being allocated on the heap as opposed to living on the stack as the individual numbers are.  
Using a static Array which is stack-allocated does not remove this difference in performance:

```julia
julia> using Statistics,StaticArrays,BenchmarkTools

julia> x = rand(Float32,10000000);StaticInput = SA[0.02,0.98];

julia> @btime quantile!($x,$StaticInput)
  490.895 ms (0 allocations: 0 bytes)
2-element SVector{2, Float64} with indices SOneTo(2):
 0.0199735164642334
 0.9800283956527709

julia> @btime quantile!($x,0.02)
  42.356 ms (0 allocations: 0 bytes)
0.0199735164642334

julia> @btime quantile!($x,0.98)
  28.156 ms (0 allocations: 0 bytes)
0.9800283956527709

```

---

<div class="post-metadata">

**Author:** ![Sukera](https://avatars.discourse-cdn.com/v4/letter/s/ce7236/32.png) [@Sukera](https://discourse.julialang.org/u/Sukera)\
**Post date:** [July 23, 2021, 4:32pm UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188/4 "2021-07-23T16:32:20Z")

</div>

Your benchmark is broken - after the first warmup, even an interpolated `x` will already be sorted because you reuse that variable all throughout the benchmark. You have to use a new `x` for each run, so use `setup=(x=rand(Float32,10000000))` and `evals=1` as additional arguments to make sure each unsorted `x` is only used once.

This doesn’t have anything to do with whether the array is allocated on the heap or not.

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [July 23, 2021, 8:21pm UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188/5 "2021-07-23T20:21:20Z")

</div>

Sorry to have brought a staring confusion by using the !inplace form. But also the mutating form showes the same behavior

```julia
x = rand(Float32, 10000000);
xx = deepcopy(x);

@time quantile(x, [0.02, 0.98])
  0.944558 seconds (4 allocations: 38.147 MiB)
2-element Vector{Float64}:
 0.019991873502731322
 0.9799974584579467

all(x .== xx)
true

@time quantile(x, 0.02)
  0.082169 seconds (3 allocations: 38.147 MiB)
0.019991873502731322

all(x .== xx)
true

@time quantile(x, 0.98)
  0.088702 seconds (3 allocations: 38.147 MiB)
0.9799974584579467

```

but also `quantil!`. See

```julia
x = rand(Float32, 10000000);

@time quantile!(x, [0.02, 0.98])
  0.947875 seconds (2 allocations: 160 bytes)
2-element Vector{Float64}:
 0.020006891489028934
 0.9800460338592529

@time quantile!(x, 0.02)
  0.021764 seconds (1 allocation: 16 bytes)
0.020006891489028934

@time quantile!(x, 0.98)
  0.021982 seconds (1 allocation: 16 bytes)
0.9800460338592529

```

and if we repeat the vector form, it’s much slower again

```julia

@time quantile!(x, [0.02, 0.98])
  0.194401 seconds (2 allocations: 160 bytes)
2-element Vector{Float64}:
 0.020006891489028934
 0.9800460338592529

```

---

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [July 23, 2021, 10:02pm UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188/6 "2021-07-23T22:02:02Z")

</div>

If we broadcast `quantile`, the times look similar to the sum of the individual times:

```julia
julia> @time quantile.(Ref(x), [0.02, 0.98])
  0.194654 seconds (9 allocations: 152.588 MiB)

# or, using BenchmarkTools:
@btime quantile.(Ref($x), [0.02, 0.98])
187.388 ms (6 allocations: 152.59 MiB)

```

---

<div class="post-metadata">

**Author:** ![genkuroki](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/genkuroki/32/18030_2.png) [@genkuroki](https://discourse.julialang.org/u/genkuroki)\
**Post date:** [July 24, 2021, 4:15am UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188/7 "2021-07-24T04:15:14Z")

</div>

**Conclusion.** The reason why `quantile!(x, [0.02, 0.98])` is much slower than `[quantile!(x, 0.02), quantile!(x, 0.98)]` is that the former executes `Statistics._quantilesort!(v, sorted, minp, maxp)` with minp \<\< maxp and the latter `Statistics._quantilesort!(v, sorted, p, p)`.

See [https://github.com/JuliaLang/Statistics.jl/blob/master/src/Statistics.jl#L949-L959](https://github.com/JuliaLang/Statistics.jl/blob/master/src/Statistics.jl#L949-L959)

**Demonstration.** ([Jupyter notebook](https://github.com/genkuroki/public/blob/main/0013/quantile!%20problem.ipynb))

```julia
using Statistics
X = rand(Float32, 10^7)

x = copy(X); y = quantile!(x, [0.02, 0.98])
x = copy(X); z = [quantile!(x, 0.02), quantile!(x, 0.98)]
@show y == z;

```

```julia
y == z = true

```

```julia
x = copy(X); @time quantile!(x, [0.02, 0.98])
x = copy(X); @time quantile!(x, [0.02, 0.98])
x = copy(X); @time quantile!(x, [0.02, 0.98]);

```

```julia
  1.171875 seconds (2 allocations: 192 bytes)
  1.175514 seconds (2 allocations: 192 bytes)
  1.174704 seconds (2 allocations: 192 bytes)

```

```julia
x = copy(X); @time [quantile!(x, 0.02), quantile!(x, 0.98)]
x = copy(X); @time [quantile!(x, 0.02), quantile!(x, 0.98)]
x = copy(X); @time [quantile!(x, 0.02), quantile!(x, 0.98)];

```

```julia
  0.179641 seconds (3 allocations: 128 bytes)
  0.196316 seconds (3 allocations: 128 bytes)
  0.178824 seconds (3 allocations: 128 bytes)

```

Yes! The former is much slower than the latter.

`@time`’s added version of [https://github.com/JuliaLang/Statistics.jl/blob/master/src/Statistics.jl#L949-L959](https://github.com/JuliaLang/Statistics.jl/blob/master/src/Statistics.jl#L949-L959)

```julia
function myquantile!(v::AbstractVector, p::Union{AbstractArray, Tuple{Vararg{Real}}};
                   sorted::Bool=false, alpha::Real=1., beta::Real=alpha)
    if !isempty(p)
        minp, maxp = extrema(p)
        print("Statistics._quantilesort!(v, sorted, minp, maxp): ")
        @time Statistics._quantilesort!(v, sorted, minp, maxp)
    end
    print("map(x->Statistics._quantile(v, x, alpha=alpha, beta=beta), p):")
    return @time map(x->Statistics._quantile(v, x, alpha=alpha, beta=beta), p)
end

myquantile!(v::AbstractVector, p::Real; sorted::Bool=false, alpha::Real=1., beta::Real=alpha) = begin
    print("Statistics._quantilesort!(v, sorted, p, p): ")
    @time Statistics._quantilesort!(v, sorted, p, p)
    print("Statistics._quantile(v, p, alpha=alpha, beta=beta):")
    @time Statistics._quantile(v, p, alpha=alpha, beta=beta)
end

println("myquantile!(x, [0.02, 0.98]):")
x = copy(X); u = myquantile!(x, [0.02, 0.98])
println("\n[myquantile!(x, 0.02), myquantile!(x, 0.98)]:")
x = copy(X); v = [myquantile!(x, 0.02), myquantile!(x, 0.98)]
println(); @show u == v;

```

Output:

```julia
myquantile!(x, [0.02, 0.98]):
Statistics._quantilesort!(v, sorted, minp, maxp): 1.208875 seconds
map(x->Statistics._quantile(v, x, alpha=alpha, beta=beta), p): 0.000001 seconds (1 allocation: 96 bytes)

[myquantile!(x, 0.02), myquantile!(x, 0.98)]:
Statistics._quantilesort!(v, sorted, p, p): 0.055211 seconds
Statistics._quantile(v, p, alpha=alpha, beta=beta): 0.000001 seconds
Statistics._quantilesort!(v, sorted, p, p): 0.098211 seconds
Statistics._quantile(v, p, alpha=alpha, beta=beta): 0.000022 seconds

u == v = true

```

As above, most of the computation time of the former is the execution time of `Statistics._quantilesort!(v, sorted, minp, maxp)` with minp \<\< maxp.

See also [https://github.com/JuliaLang/julia/blob/361b3a48e2ce05a05f7e3c01903d18bf114432e3/base/sort.jl#L628-L652](https://github.com/JuliaLang/julia/blob/361b3a48e2ce05a05f7e3c01903d18bf114432e3/base/sort.jl#L628-L652)

**Postscript.** (case of [0.499, 0.501])

```julia
x = copy(X); @time y = quantile!(x, [0.499, 0.501])
x = copy(X); @time y = quantile!(x, [0.499, 0.501])
x = copy(X); @time y = quantile!(x, [0.499, 0.501]);

```

```julia
  0.133984 seconds (2 allocations: 192 bytes)
  0.138859 seconds (2 allocations: 192 bytes)
  0.161751 seconds (2 allocations: 192 bytes)

```

```julia
x = copy(X); @time z = [quantile!(x, 0.499), quantile!(x, 0.501)]
x = copy(X); @time z = [quantile!(x, 0.499), quantile!(x, 0.501)]
x = copy(X); @time z = [quantile!(x, 0.499), quantile!(x, 0.501)];

```

```julia
  0.235203 seconds (3 allocations: 128 bytes)
  0.225585 seconds (3 allocations: 128 bytes)
  0.200007 seconds (3 allocations: 128 bytes)

```

In this case, the former is faster than the latter.

---

<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:** [August 5, 2021, 9:58pm UTC](https://discourse.julialang.org/t/quantile-much-faster-when-called-twice-with-a-scalar-than-once-with-a-vector/65188/8 "2021-08-05T21:58:33Z")

</div>

[EmpiricalCDFs.jl](https://github.com/jlapeyre/EmpiricalCDFs.jl) uses a different approach that may be more efficient, depending on your use case

```julia
julia> using EmpiricalCDFs

julia> X = rand(Float32, 10^7);

julia> cdf = EmpiricalCDF(X);

julia> sort!(cdf);

julia> const icdf2 = finv(cdf); # get the quantile function

julia> @btime icdf2.([.02, .98])
  42.028 ns (2 allocations: 192 bytes)
2-element Vector{Float32}:
 0.02000904
 0.9799932

```