stat: add support for LinInterp CumulantKind for Quantile
MODIFY
stat/stat.go
stat/stat.go
@@ -20,6 +20,8 @@ type CumulantKind int
const (
// Empirical treats the distribution as the actual empirical distribution.
Empirical CumulantKind = 1
+ // LinInterp linearly interpolates the empirical distribution between sample values, with a flat extrapolation.
+ LinInterp CumulantKind = 4
)
// bhattacharyyaCoeff computes the Bhattacharyya Coefficient for probability distributions given by:
@@ -1028,6 +1030,7 @@ func MomentAbout(moment float64, x []float64, mean float64, weights []float64) f
// CumulantKind behaviors:
// - Empirical: Returns the lowest value q for which q is greater than or equal
// to the fraction p of samples
+// - LinInterp: Returns the linearly interpolated value
func Quantile(p float64, c CumulantKind, x, weights []float64) float64 {
if !(p >= 0 && p <= 1) {
panic("stat: percentile out of bounds")
@@ -1052,6 +1055,8 @@ func Quantile(p float64, c CumulantKind, x, weights []float64) float64 {
switch c {
case Empirical:
return empiricalQuantile(p, x, weights, sumWeights)
+ case LinInterp:
+ return linInterpQuantile(p, x, weights, sumWeights)
default:
panic("stat: bad cumulant kind")
}
@@ -1073,6 +1078,29 @@ func empiricalQuantile(p float64, x, weights []float64, sumWeights float64) floa
panic("impossible")
}
+func linInterpQuantile(p float64, x, weights []float64, sumWeights float64) float64 {
+ var cumsum float64
+ fidx := p * sumWeights
+ for i := range x {
+ if weights == nil {
+ cumsum++
+ } else {
+ cumsum += weights[i]
+ }
+ if cumsum >= fidx {
+ if i == 0 {
+ return x[0]
+ }
+ t := cumsum - fidx
+ if weights != nil {
+ t /= weights[i]
+ }
+ return t*x[i-1] + (1-t)*x[i]
+ }
+ }
+ panic("impossible")
+}
+
// Skew computes the skewness of the sample data.
// If weights is nil then all of the weights are 1. If weights is not nil, then
// len(x) must equal len(weights).

MODIFY
stat/stat_test.go
stat/stat_test.go
@@ -1423,7 +1423,10 @@ func TestCDF(t *testing.T) {
}
func TestQuantile(t *testing.T) {
- cumulantKinds := []CumulantKind{Empirical}
+ cumulantKinds := []CumulantKind{
+ Empirical,
+ LinInterp,
+ }
for i, test := range []struct {
p []float64
x []float64
@@ -1434,18 +1437,36 @@ func TestQuantile(t *testing.T) {
p: []float64{0, 0.05, 0.1, 0.15, 0.45, 0.5, 0.55, 0.85, 0.9, 0.95, 1},
x: []float64{1, 2, 3, 4, 5, 6, 7, 8, 9, 10},
w: nil,
- ans: [][]float64{{1, 1, 1, 2, 5, 5, 6, 9, 9, 10, 10}},
+ ans: [][]float64{
+ {1, 1, 1, 2, 5, 5, 6, 9, 9, 10, 10},
+ {1, 1, 1, 1.5, 4.5, 5, 5.5, 8.5, 9, 9.5, 10},
+ },
},
{
p: []float64{0, 0.05, 0.1, 0.15, 0.45, 0.5, 0.55, 0.85, 0.9, 0.95, 1},
x: []float64{1, 2, 3, 4, 5, 6, 7, 8, 9, 10},
w: []float64{3, 3, 3, 3, 3, 3, 3, 3, 3, 3},
- ans: [][]float64{{1, 1, 1, 2, 5, 5, 6, 9, 9, 10, 10}},
+ ans: [][]float64{
+ {1, 1, 1, 2, 5, 5, 6, 9, 9, 10, 10},
+ {1, 1, 1, 1.5, 4.5, 5, 5.5, 8.5, 9, 9.5, 10},
+ },
+ },
+ {
+ p: []float64{0, 0.05, 0.1, 0.15, 0.45, 0.5, 0.55, 0.85, 0.9, 0.95, 1},
+ x: []float64{1, 2, 3, 4, 5, 6, 7, 8, 9, 10},
+ w: []float64{0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0},
+ ans: [][]float64{
+ {1, 2, 3, 4, 7, 7, 8, 10, 10, 10, 10},
+ {1, 1.875, 2.833333333333333, 3.5625, 6.535714285714286, 6.928571428571429, 7.281250000000001, 9.175, 9.45, 9.725, 10},
+ },
},
{
p: []float64{0.5},
x: []float64{1, 2, 3, 4, 5, 6, 7, 8, math.NaN(), 10},
- ans: [][]float64{{math.NaN()}},
+ ans: [][]float64{
+ {math.NaN()},
+ {math.NaN()},
+ },
},
} {
copyX := make([]float64, len(test.x))
@@ -1470,54 +1491,56 @@ func TestQuantile(t *testing.T) {
}
}
}
- // panic cases
+}
+
+func TestQuantileInvalidInput(t *testing.T) {
+ cumulantKinds := []CumulantKind{
+ Empirical,
+ LinInterp,
+ }
for _, test := range []struct {
name string
p float64
- c CumulantKind
x []float64
w []float64
}{
{
name: "p < 0",
- c: Empirical,
p: -1,
},
{
name: "p > 1",
- c: Empirical,
p: 2,
},
{
name: "p is NaN",
- c: Empirical,
p: math.NaN(),
},
{
name: "len(x) != len(weights)",
- c: Empirical,
p: .5,
x: make([]float64, 4),
w: make([]float64, 2),
},
{
name: "x not sorted",
- c: Empirical,
p: .5,
x: []float64{3, 2, 1},
},
- {
- name: "CumulantKind is unknown",
- c: CumulantKind(1000),
- p: .5,
- x: []float64{1, 2, 3},
- },
} {
- if !panics(func() { Quantile(test.p, test.c, test.x, test.w) }) {
+ for _, kind := range cumulantKinds {
+ if !panics(func() { Quantile(test.p, kind, test.x, test.w) }) {
t.Errorf("Quantile did not panic when %s", test.name)
}
}
}
+}
+
+func TestQuantileInvalidCumulantKind(t *testing.T) {
+ if !panics(func() { Quantile(0.5, CumulantKind(1000), []float64{1, 2, 3}, nil) }) {
+ t.Errorf("Quantile did not panic when CumulantKind is unknown")
+ }
+}
func ExampleStdDev() {
x := []float64{8, 2, -9, 15, 4}

stat: add support for LinInterp CumulantKind for Quantile