Merge pull request #16572 from JuliaLang/sb/quantile-acc

Fixes numerical accuracy issues in quantile.

base/statistics.jl

@@ -577,19 +577,21 @@ end
 @inline function _quantile(v::AbstractVector, p::Real)
  T = float(eltype(v))
  isnan(p) && return T(NaN)
+ 0 <= p <= 1 || throw(ArgumentError("input probability out of [0,1] range"))
 
  lv = length(v)
- index = 1 + (lv-1)*p
- 1 <= index <= lv || error("input probability out of [0,1] range")
+ f0 = (lv-1)*p # 0-based interpolated index
+ t0 = trunc(f0)
+ h = f0 - t0
 
- indlo = floor(index)
- i = trunc(Int,indlo)
+ i = trunc(Int,t0) + 1
 
- if index == indlo
+ if h == 0
  return T(v[i])
  else
- h = T(index - indlo)
- return (1-h)*T(v[i]) + h*T(v[i+1])
+ a = T(v[i])
+ b = T(v[i+1])
+ return a + h*(b-a)
  end
 end
 

test/statistics.jl

@@ -324,6 +324,11 @@ end
 @test quantile(0.0:100.0, 0.0:0.1:1.0, sorted=true) == collect(0.0:10.0:100.0)
 @test quantile(100f0:-1f0:0.0, 0.0:0.1:1.0) == collect(0f0:10f0:100f0)
 
+@test quantile([0,1],1e-18) == 1e-18
+
+# StatsBase issue 164
+y = [0.40003674665581906,0.4085630862624367,0.41662034698690303,0.41662034698690303,0.42189053966652057,0.42189053966652057,0.42553514344518345,0.43985732442991354]
+@test issorted(quantile(y, linspace(0.01, 0.99, 17)))
 
 # test invalid hist nbins argument (#9999)
 @test_throws ArgumentError hist(Int[], -1)