Fix median for Binomial distribution (#1924)

* Fix `median` for `Binomial` Fix #1923 * Add tests for `median(x::Binomial)` * Remove specialization of `median` to `Binomial` Let `median(x::Binomial)` fall back to `quantile(x,1//2)` * Update R reference * Simplify test --------- Co-authored-by: David Widmann <[email protected]>
JuliaStats · Dec 14, 2024 · 2d28655 · 2d28655
1 parent 3de6038
commit 2d28655
Show file tree

Hide file tree

Showing 4 changed files with 9 additions and 4 deletions.
diff --git a/src/univariate/discrete/binomial.jl b/src/univariate/discrete/binomial.jl
@@ -73,8 +73,6 @@ function mode(d::Binomial{T}) where T<:Real
 end
 modes(d::Binomial) = Int[mode(d)]
 
-median(d::Binomial) = round(Int,mean(d))
-
 function skewness(d::Binomial)
     n, p1 = params(d)
     p0 = 1 - p1

diff --git a/test/ref/discrete/binomial.R b/test/ref/discrete/binomial.R
@@ -18,7 +18,7 @@ Binomial <- R6Class("Binomial",
                  failprob=q,
                  ntrials=n,
                  mean=n * p,
-                 median=round(n * p),
+                 median=qbinom(0.5, self$n, self$p),
                  var=n * p * q,
                  skewness=(q - p) / sqrt(n*p*q),
                  kurtosis=(1 - 6*p*q) / (n*p*q))

diff --git a/test/ref/discrete_test.ref.json b/test/ref/discrete_test.ref.json
@@ -285,7 +285,7 @@
     "failprob": 0.5,
     "ntrials": 3,
     "mean": 1.5,
-    "median": 2,
+    "median": 1,
     "var": 0.75,
     "skewness": 0,
     "kurtosis": -0.666666666666667

diff --git a/test/univariate/discrete/binomial.jl b/test/univariate/discrete/binomial.jl
@@ -26,6 +26,13 @@ end
 @test Distributions.expectation(identity, Binomial(6)) ≈ 3.0
 @test Distributions.expectation(x -> -x, Binomial(10, 0.2)) ≈ -2.0
 
+# Test median
+@test median(Binomial(5,3//10)) == 1
+@test median(Binomial(25,3//10)) == 7
+@test median(Binomial(45,3//10)) == 13
+@test median(Binomial(65,3//10)) == 19
+@test median(Binomial(85,3//10)) == 25
+
 # Test mode
 @test Distributions.mode(Binomial(100, 0.4)) == 40
 @test Distributions.mode(Binomial(1, 0.51)) == 1