diff --git a/src/univariate/discrete/binomial.jl b/src/univariate/discrete/binomial.jl
index 74c90ee18..5390dbce1 100644
--- a/src/univariate/discrete/binomial.jl
+++ b/src/univariate/discrete/binomial.jl
@@ -67,11 +67,36 @@ params(d::Binomial) = (d.n, d.p)
 
 mean(d::Binomial) = ntrials(d) * succprob(d)
 var(d::Binomial) = ntrials(d) * succprob(d) * failprob(d)
-function mode(d::Binomial{T}) where T<:Real
+
+function mode(d::Binomial)
+    (n, p) = params(d)
+    v = (n + 1) * p
+    quasi_mode = floor(Int, v)
+    if quasi_mode == v && p > 0
+        if p == 1
+            n
+        else
+            quasi_mode-1
+        end
+    else
+        quasi_mode
+    end
+end
+
+function modes(d::Binomial)
     (n, p) = params(d)
-    n > 0 ? floor(Int, (n + 1) * d.p) : zero(T)
+    v = (n + 1) * p
+    quasi_mode = floor(Int, v)
+    if quasi_mode == v
+        if p == 1
+            Int[n]
+        else
+            Int[quasi_mode-1, quasi_mode]
+        end
+    else
+        Int[quasi_mode]
+    end
 end
-modes(d::Binomial) = Int[mode(d)]
 
 function median(dist::Binomial)
     # The median is floor(Int, mean) or ceil(Int, mean)
diff --git a/test/univariate/discrete/binomial.jl b/test/univariate/discrete/binomial.jl
index 55b336214..3b49f30b6 100644
--- a/test/univariate/discrete/binomial.jl
+++ b/test/univariate/discrete/binomial.jl
@@ -35,9 +35,20 @@ end
 @test all(median(Binomial(7, p)) == quantile(Binomial(7, p), 1//2) for p in 0:0.1:1)
     
 # Test mode
-@test Distributions.mode(Binomial(100, 0.4)) == 40
-@test Distributions.mode(Binomial(1, 0.51)) == 1
-@test Distributions.mode(Binomial(1, 0.49)) == 0
+@test mode(Binomial(100, 0.4)) == 40
+@test mode(Binomial(1, 0.51)) == 1
+@test mode(Binomial(1, 0.49)) == 0
+@test mode(Binomial(4, 2//3)) == 3
+@test mode(Binomial(6, 2//7)) == 1
+@test mode(Binomial(7, 1//8)) == 0
+
+@test modes(Binomial(4, 2//3)) == [3]
+@test modes(Binomial(5, 3//4)) == [4]
+@test modes(Binomial(3, 2//4)) == [1, 2]
+@test modes(Binomial(4, 2//5)) == [1, 2]
+@test modes(Binomial(6, 2//7)) == [1, 2]
+@test modes(Binomial(6, 3//7)) == [2, 3]
+@test modes(Binomial(7, 1//8)) == [0, 1]
 
 @test isplatykurtic(Bernoulli(0.5))
 @test ismesokurtic(Normal(0.0, 1.0))