quantile.rb

#!/usr/bin/env ruby

### unsorted data series
Y = [
    0.049632, 0.056326, 0.052863, 0.054129, 0.019590,
    0.010500, 0.038988, 0.276094, 0.028402, 0.030384,
    0.038914, 0.029160, 0.066111, 0.015465, 0.076519,
    0.048751, 0.043617, 0.005298, 0.005722, 0.049853,
    0.059949, 0.009066, 0.442716, 0.012065, 0.029591,
    0.053335, 0.048495, 0.059594, 0.061983, 0.014104,
    0.358289, 0.025765, 0.054330, 0.073260, 0.061888,
    0.066990, 0.059203, 0.076802, 0.015893, 0.026751,
    0.061616, 0.072682, 0.035666, 0.066835, 0.054902,
    0.072753, 0.190644, 0.036291, 0.065030, 0.076660
]

    ### Lagrange polynomial interpolation
    ### ref: https://en.wikipedia.org/wiki/Lagrange_polynomial

    ### Lagrange polynomial interpolation: helper function
    def lx(f, x, i)
        ret = 1.0
        f.count.times { |k|
            next if k == i
            ret *= (x - f[k][0])
        }
        f.count.times { |k|
            next if k == i
            ret /= (f[i][0] - f[k][0])
        }
        return ret
    end

    ### Lagrange polynomial interpolation: main function
    def Lx(f, x)
        ret = 0.0
        f.count.times { |i|
            ret += (f[i][1] * lx(f, x, i))
        }
        return ret
    end

    def Quantile(f, a)
        ### select 2 adjacent "x->y" pairs near quantile value
        z = [ f.select { |xy| xy[0] <= a }.last, f.select { |xy| xy[0] > a }.first ]

    ##  return Lx(z, a)

        ### data series with 2 values doesn't need entire Lagrange polynomial
        ### here is plain formula:
        ###   q = (a*(y1-y0)+x1*y0+x0*y1)/(x1-x0)
        ret = 0.0
        ret += (z[1][0] * z[0][1] - z[0][0] * z[1][1])
        ret += a * (z[1][1] - z[0][1])
        ret /= (z[1][0] - z[0][0])

        return ret
    end

### in-place sort
Y.sort!

### create array of "x->y" pairs (via array with two element)
###   x = normalized index in [0;1]
###   y = sorted data series
F = []
Y.count.times { |z| F << [ z / (Y.count - 1.0), Y[z] ] }

### kind of "pretty output"
# F.count.times { |i|
#   printf "%-10f :: %f\n", F[i][0], F[i][1]
# }
# puts

### calculate several quantiles
[0.05, 0.10, 0.15, 0.20, 0.25, 0.50, 0.75, 0.80, 0.85, 0.90, 0.95].each { |q|
    z = Quantile(F, q)
    printf "%-10f => %f (%d)\n", q, z, F.count { |xy| xy[1] < z }
}

# $ ./quantile.rb
# 0.000000   :: 0.005298
# 0.020408   :: 0.005722
# 0.040816   :: 0.009066
# 0.061224   :: 0.010500
# 0.081633   :: 0.012065
# 0.102041   :: 0.014104
# 0.122449   :: 0.015465
# 0.142857   :: 0.015893
# 0.163265   :: 0.019590
# 0.183673   :: 0.025765
# 0.204082   :: 0.026751
# 0.224490   :: 0.028402
# 0.244898   :: 0.029160
# 0.265306   :: 0.029591
# 0.285714   :: 0.030384
# 0.306122   :: 0.035666
# 0.326531   :: 0.036291
# 0.346939   :: 0.038914
# 0.367347   :: 0.038988
# 0.387755   :: 0.043617
# 0.408163   :: 0.048495
# 0.428571   :: 0.048751
# 0.448980   :: 0.049632
# 0.469388   :: 0.049853
# 0.489796   :: 0.052863
# 0.510204   :: 0.053335
# 0.530612   :: 0.054129
# 0.551020   :: 0.054330
# 0.571429   :: 0.054902
# 0.591837   :: 0.056326
# 0.612245   :: 0.059203
# 0.632653   :: 0.059594
# 0.653061   :: 0.059949
# 0.673469   :: 0.061616
# 0.693878   :: 0.061888
# 0.714286   :: 0.061983
# 0.734694   :: 0.065030
# 0.755102   :: 0.066111
# 0.775510   :: 0.066835
# 0.795918   :: 0.066990
# 0.816327   :: 0.072682
# 0.836735   :: 0.072753
# 0.857143   :: 0.073260
# 0.877551   :: 0.076519
# 0.897959   :: 0.076660
# 0.918367   :: 0.076802
# 0.938776   :: 0.190644
# 0.959184   :: 0.276094
# 0.979592   :: 0.358289
# 1.000000   :: 0.442716
#
# 0.050000   => 0.009711 (3)
# 0.100000   => 0.013900 (5)
# 0.150000   => 0.017187 (8)
# 0.200000   => 0.026554 (10)
# 0.250000   => 0.029268 (13)
# 0.500000   => 0.053099 (25)
# 0.750000   => 0.065841 (37)
# 0.800000   => 0.068128 (40)
# 0.850000   => 0.073083 (42)
# 0.900000   => 0.076674 (45)
# 0.950000   => 0.237642 (47)