@@ -89,7 +89,7 @@ namespace MathNet.Numerics.Distributions /// <returns>false <paramref name="p"/> is not in the interval [0.0,1.0] or <paramref name="n"/> is negative, true otherwise.</exception> private static bool IsValidParameterSet(double p, int n) { - if(p < 0.0 || p > 1.0) + if(p < 0.0 || p > 1.0 || Double.IsNaN(p)) { return false; } @@ -206,11 +206,16 @@ namespace MathNet.Numerics.Distributions { get { + if (_p == 0.0 || _p == 1.0) + { + return 0.0; + } + double E = 0.0; - for(int i = 0; i < _n; i++) + for(int i = 0; i <= _n; i++) { double p = Probability(i); - E += p * Math.Log(p); + E -= p * Math.Log(p); } return E; } @@ -263,7 +268,19 @@ namespace MathNet.Numerics.Distributions /// </summary> public int Mode { - get { return (int) Math.Floor((_n + 1) * _p); } + get + { + if (_p == 1.0) + { + return _n; + } + else if (_p == 0.0) + { + return 0; + } + + return (int) Math.Floor((_n + 1) * _p); + } } /// <summary> @@ -289,23 +306,59 @@ namespace MathNet.Numerics.Distributions return 0.0; } + if (_p == 0.0 && val == 0) + { + return 1.0; + } + else if (_p == 0.0) + { + return 0.0; + } + + if (_p == 1.0 && val == _n) + { + return 1.0; + } + else if (_p == 1.0) + { + return 0.0; + } + return SpecialFunctions.Binomial(_n, val) * Math.Pow(_p, val) * Math.Pow(1.0 - _p, _n - val); } /// <summary> - /// Computes the probability of a specific value. + /// Computes the log probability of a specific value. /// </summary> public double ProbabilityLn(int val) { if (val < 0) { - return 0.0; + return Double.NegativeInfinity; } if (val > _n) { + return Double.NegativeInfinity; + } + + if (_p == 0.0 && val == 0) + { return 0.0; } + else if (_p == 0.0) + { + return Double.NegativeInfinity; + } + + if (_p == 1.0 && val == _n) + { + return 0.0; + } + else if (_p == 1.0) + { + return Double.NegativeInfinity; + } return SpecialFunctions.BinomialLn(_n, val) + val * Math.Log(_p) + (_n - val) * Math.Log(1.0 - _p); } src/Numerics/Distributions/Discrete/Binomial.cs @@ -52,7 +52,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests dists[6] = new DiscreteUniform(1, 10); dists[7] = new LogNormal(1.0, 1.0); dists[8] = new Binomial(0.7, 10); - dists[9] = new Categorical(0.7); + dists[9] = new Categorical(new double[] { 0.7, 0.3 }); } [Test] src/UnitTests/DistributionTests/CommonDistributionTests.cs @@ -92,13 +92,13 @@ namespace MathNet.Numerics.UnitTests.DistributionTests } [Test] - [Row(0.0, 4)] - [Row(0.3, 3)] - [Row(1.0, 2)] - public void ValidateEntropy(double p, int n) + [Row(0.0, 4, 0.0)] + [Row(0.3, 3, 1.1404671643037712668976423399228972051669206536461)] + [Row(1.0, 2, 0.0)] + public void ValidateEntropy(double p, int n, double e) { var b = new Binomial(p,n); - AssertHelpers.AlmostEqual(n * (-(1.0 - p) * Math.Log(1.0 - p) - p * Math.Log(p)), b.Entropy, 14); + AssertHelpers.AlmostEqual(e, b.Entropy, 14); } [Test] @@ -112,13 +112,13 @@ namespace MathNet.Numerics.UnitTests.DistributionTests } [Test] - [Row(0.0, 4, 0.0)] - [Row(0.3, 3, 1.0)] - [Row(1.0, 2, 1.0)] - public void ValidateMode(double p, int n, double m) + [Row(0.0, 4, 0)] + [Row(0.3, 3, 1)] + [Row(1.0, 2, 2)] + public void ValidateMode(double p, int n, int m) { var b = new Binomial(p,n); - AssertEx.AreEqual<double>(m, b.Mode); + AssertEx.AreEqual<int>(m, b.Mode); } [Test] @@ -138,7 +138,6 @@ namespace MathNet.Numerics.UnitTests.DistributionTests [Test] [Row(0.000000, 1, 0, 1.0)] [Row(0.000000, 1, 1, 0.0)] - [Row(0.000000, 1, 1, 0.0)] [Row(0.000000, 3, 0, 1.0)] [Row(0.000000, 3, 1, 0.0)] [Row(0.000000, 3, 3, 0.0)] @@ -147,7 +146,6 @@ namespace MathNet.Numerics.UnitTests.DistributionTests [Row(0.000000, 10, 10, 0.0)] [Row(0.300000, 1, 0, 0.69999999999999995559107901499373838305473327636719)] [Row(0.300000, 1, 1, 0.2999999999999999888977697537484345957636833190918)] - [Row(0.300000, 1, 1, 0.2999999999999999888977697537484345957636833190918)] [Row(0.300000, 3, 0, 0.34299999999999993471888615204079956461021032657166)] [Row(0.300000, 3, 1, 0.44099999999999992772448109690231306411849135972008)] [Row(0.300000, 3, 3, 0.026999999999999997002397833512077451789759292859569)] @@ -156,7 +154,6 @@ namespace MathNet.Numerics.UnitTests.DistributionTests [Row(0.300000, 10, 10, 0.0000059048999999999978147480206303047454017251032868501)] [Row(1.000000, 1, 0, 0.0)] [Row(1.000000, 1, 1, 1.0)] - [Row(1.000000, 1, 1, 1.0)] [Row(1.000000, 3, 0, 0.0)] [Row(1.000000, 3, 1, 0.0)] [Row(1.000000, 3, 3, 1.0)] @@ -166,41 +163,38 @@ namespace MathNet.Numerics.UnitTests.DistributionTests public void ValidateProbability(double p, int n, int x, double d) { var b = new Binomial(p,n); - AssertEx.AreEqual(d, b.Probability(x)); + AssertHelpers.AlmostEqual(d, b.Probability(x), 14); } [Test] [Row(0.000000, 1, 0, 0.0)] - [Row(0.000000, 1, 1, -inf)] - [Row(0.000000, 1, 1, -inf)] + [Row(0.000000, 1, 1, Double.NegativeInfinity)] [Row(0.000000, 3, 0, 0.0)] - [Row(0.000000, 3, 1, -inf)] - [Row(0.000000, 3, 3, -inf)] + [Row(0.000000, 3, 1, Double.NegativeInfinity)] + [Row(0.000000, 3, 3, Double.NegativeInfinity)] [Row(0.000000, 10, 0, 0.0)] - [Row(0.000000, 10, 1, -inf)] - [Row(0.000000, 10, 10, -inf)] + [Row(0.000000, 10, 1, Double.NegativeInfinity)] + [Row(0.000000, 10, 10, Double.NegativeInfinity)] [Row(0.300000, 1, 0, -0.3566749439387324423539544041072745145718090708995)] [Row(0.300000, 1, 1, -1.2039728043259360296301803719337238685164245381839)] - [Row(0.300000, 1, 1, -1.2039728043259360296301803719337238685164245381839)] [Row(0.300000, 3, 0, -1.0700248318161973270618632123218235437154272126985)] [Row(0.300000, 3, 1, -0.81871040353529122294284394322574719301255212216016)] [Row(0.300000, 3, 3, -3.6119184129778080888905411158011716055492736145517)] [Row(0.300000, 10, 0, -3.566749439387324423539544041072745145718090708995)] [Row(0.300000, 10, 1, -2.1114622067804823267977785542148302920616046876506)] [Row(0.300000, 10, 10, -12.039728043259360296301803719337238685164245381839)] - [Row(1.000000, 1, 0, -inf)] - [Row(1.000000, 1, 1, 0.0)] + [Row(1.000000, 1, 0, Double.NegativeInfinity)] [Row(1.000000, 1, 1, 0.0)] - [Row(1.000000, 3, 0, -inf)] - [Row(1.000000, 3, 1, -inf)] + [Row(1.000000, 3, 0, Double.NegativeInfinity)] + [Row(1.000000, 3, 1, Double.NegativeInfinity)] [Row(1.000000, 3, 3, 0.0)] - [Row(1.000000, 10, 0, -inf)] - [Row(1.000000, 10, 1, -inf)] + [Row(1.000000, 10, 0, Double.NegativeInfinity)] + [Row(1.000000, 10, 1, Double.NegativeInfinity)] [Row(1.000000, 10, 10, 0.0)] public void ValidateProbabilityLn(double p, int n, int x, double dln) { var b = new Binomial(p,n); - AssertEx.AreEqual(dln, b.ProbabilityLn(x)); + AssertHelpers.AlmostEqual(dln, b.ProbabilityLn(x), 14); } [Test] src/UnitTests/DistributionTests/Discrete/BinomialTests.cs @@ -54,7 +54,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests [Test] public void CanCreateCategorical() { - var m = new Categorical(largeP, 4); + var m = new Categorical(largeP); AssertEx.AreEqual<double[]>(largeP, m.P); } @@ -62,14 +62,14 @@ namespace MathNet.Numerics.UnitTests.DistributionTests [ExpectedException(typeof(ArgumentOutOfRangeException))] public void CategoricalCreateFailsWithNegativeRatios() { - var m = new Categorical(badP, 4); + var m = new Categorical(badP); } [Test] [ExpectedException(typeof(ArgumentOutOfRangeException))] public void CategoricalCreateFailsWithAllZeroRatios() { - var m = new Categorical(badP2, 4); + var m = new Categorical(badP2); } [Test] @@ -82,7 +82,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests [Test] public void CanSetProbability() { - var b = new Categorical(largeP, 4); + var b = new Categorical(largeP); b.P = smallP; } @@ -90,27 +90,27 @@ namespace MathNet.Numerics.UnitTests.DistributionTests [ExpectedException(typeof(ArgumentOutOfRangeException))] public void SetProbabilityFails() { - var b = new Categorical(largeP, 4); + var b = new Categorical(largeP); b.P = badP; } [Test] public void CanSampleStatic() { - var d = Categorical.Sample(new Random(), largeP, 4); + var d = Categorical.Sample(new Random(), largeP); } [Test] [ExpectedException(typeof(ArgumentOutOfRangeException))] public void FailSampleStatic() { - var d = Categorical.Sample(new Random(), badP, 4); + var d = Categorical.Sample(new Random(), badP); } [Test] public void CanSample() { - var n = new Categorical(largeP, 4); + var n = new Categorical(largeP); var d = n.Sample(); } } src/UnitTests/DistributionTests/Discrete/CategoricalTests.cs @@ -76,7 +76,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests public void ValidateToString() { var b = new Multinomial(smallP, 4); - AssertEx.AreEqual<string>("Multinomial(Dimension = 3)", b.ToString()); + AssertEx.AreEqual<string>("Multinomial(Dimension = 3, Number of Trails = 4)", b.ToString()); } [Test] src/UnitTests/DistributionTests/Multivariate/MultinomialTests.cs 2624f484212c54ba80dceb7c777a12f27f3d97a2 jvangael jurgen.vangael@gmail.com http://github.com/cuda/mathnet-numerics/commit/96dbecc573f040e68e09e15adf55e650b7bcbccc 96dbecc573f040e68e09e15adf55e650b7bcbccc 2009-11-01T14:37:23-08:00 2009-11-01T14:20:44-08:00 Fixed bugs in distributions. b7ac26f97a4206bfd7483ef8f2bb2e5479fe8686 jvangael jurgen.vangael@gmail.com