Browse Source

Fixed bugs in distributions.

la-knuth
Jurgen Van Gael 17 years ago
parent
commit
72e771abbb
  1. 65
      src/Numerics/Distributions/Discrete/Binomial.cs
  2. 2
      src/UnitTests/DistributionTests/CommonDistributionTests.cs
  3. 50
      src/UnitTests/DistributionTests/Discrete/BinomialTests.cs
  4. 16
      src/UnitTests/DistributionTests/Discrete/CategoricalTests.cs
  5. 2
      src/UnitTests/DistributionTests/Multivariate/MultinomialTests.cs

65
src/Numerics/Distributions/Discrete/Binomial.cs

@ -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> /// <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) 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; return false;
} }
@ -206,11 +206,16 @@ namespace MathNet.Numerics.Distributions
{ {
get get
{ {
if (_p == 0.0 || _p == 1.0)
{
return 0.0;
}
double E = 0.0; double E = 0.0;
for(int i = 0; i < _n; i++) for(int i = 0; i <= _n; i++)
{ {
double p = Probability(i); double p = Probability(i);
E += p * Math.Log(p); E -= p * Math.Log(p);
} }
return E; return E;
} }
@ -263,7 +268,19 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public int Mode 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> /// <summary>
@ -289,23 +306,59 @@ namespace MathNet.Numerics.Distributions
return 0.0; 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); return SpecialFunctions.Binomial(_n, val) * Math.Pow(_p, val) * Math.Pow(1.0 - _p, _n - val);
} }
/// <summary> /// <summary>
/// Computes the probability of a specific value. /// Computes the log probability of a specific value.
/// </summary> /// </summary>
public double ProbabilityLn(int val) public double ProbabilityLn(int val)
{ {
if (val < 0) if (val < 0)
{ {
return 0.0; return Double.NegativeInfinity;
} }
if (val > _n) if (val > _n)
{
return Double.NegativeInfinity;
}
if (_p == 0.0 && val == 0)
{ {
return 0.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); return SpecialFunctions.BinomialLn(_n, val) + val * Math.Log(_p) + (_n - val) * Math.Log(1.0 - _p);
} }

2
src/UnitTests/DistributionTests/CommonDistributionTests.cs

@ -52,7 +52,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
dists[6] = new DiscreteUniform(1, 10); dists[6] = new DiscreteUniform(1, 10);
dists[7] = new LogNormal(1.0, 1.0); dists[7] = new LogNormal(1.0, 1.0);
dists[8] = new Binomial(0.7, 10); dists[8] = new Binomial(0.7, 10);
dists[9] = new Categorical(0.7); dists[9] = new Categorical(new double[] { 0.7, 0.3 });
} }
[Test] [Test]

50
src/UnitTests/DistributionTests/Discrete/BinomialTests.cs

@ -92,13 +92,13 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
} }
[Test] [Test]
[Row(0.0, 4)] [Row(0.0, 4, 0.0)]
[Row(0.3, 3)] [Row(0.3, 3, 1.1404671643037712668976423399228972051669206536461)]
[Row(1.0, 2)] [Row(1.0, 2, 0.0)]
public void ValidateEntropy(double p, int n) public void ValidateEntropy(double p, int n, double e)
{ {
var b = new Binomial(p,n); 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] [Test]
@ -112,13 +112,13 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
} }
[Test] [Test]
[Row(0.0, 4, 0.0)] [Row(0.0, 4, 0)]
[Row(0.3, 3, 1.0)] [Row(0.3, 3, 1)]
[Row(1.0, 2, 1.0)] [Row(1.0, 2, 2)]
public void ValidateMode(double p, int n, double m) public void ValidateMode(double p, int n, int m)
{ {
var b = new Binomial(p,n); var b = new Binomial(p,n);
AssertEx.AreEqual<double>(m, b.Mode); AssertEx.AreEqual<int>(m, b.Mode);
} }
[Test] [Test]
@ -138,7 +138,6 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
[Test] [Test]
[Row(0.000000, 1, 0, 1.0)] [Row(0.000000, 1, 0, 1.0)]
[Row(0.000000, 1, 1, 0.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, 0, 1.0)]
[Row(0.000000, 3, 1, 0.0)] [Row(0.000000, 3, 1, 0.0)]
[Row(0.000000, 3, 3, 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.000000, 10, 10, 0.0)]
[Row(0.300000, 1, 0, 0.69999999999999995559107901499373838305473327636719)] [Row(0.300000, 1, 0, 0.69999999999999995559107901499373838305473327636719)]
[Row(0.300000, 1, 1, 0.2999999999999999888977697537484345957636833190918)] [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, 0, 0.34299999999999993471888615204079956461021032657166)]
[Row(0.300000, 3, 1, 0.44099999999999992772448109690231306411849135972008)] [Row(0.300000, 3, 1, 0.44099999999999992772448109690231306411849135972008)]
[Row(0.300000, 3, 3, 0.026999999999999997002397833512077451789759292859569)] [Row(0.300000, 3, 3, 0.026999999999999997002397833512077451789759292859569)]
@ -156,7 +154,6 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
[Row(0.300000, 10, 10, 0.0000059048999999999978147480206303047454017251032868501)] [Row(0.300000, 10, 10, 0.0000059048999999999978147480206303047454017251032868501)]
[Row(1.000000, 1, 0, 0.0)] [Row(1.000000, 1, 0, 0.0)]
[Row(1.000000, 1, 1, 1.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, 0, 0.0)]
[Row(1.000000, 3, 1, 0.0)] [Row(1.000000, 3, 1, 0.0)]
[Row(1.000000, 3, 3, 1.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) public void ValidateProbability(double p, int n, int x, double d)
{ {
var b = new Binomial(p,n); var b = new Binomial(p,n);
AssertEx.AreEqual(d, b.Probability(x)); AssertHelpers.AlmostEqual(d, b.Probability(x), 14);
} }
[Test] [Test]
[Row(0.000000, 1, 0, 0.0)] [Row(0.000000, 1, 0, 0.0)]
[Row(0.000000, 1, 1, -inf)] [Row(0.000000, 1, 1, Double.NegativeInfinity)]
[Row(0.000000, 1, 1, -inf)]
[Row(0.000000, 3, 0, 0.0)] [Row(0.000000, 3, 0, 0.0)]
[Row(0.000000, 3, 1, -inf)] [Row(0.000000, 3, 1, Double.NegativeInfinity)]
[Row(0.000000, 3, 3, -inf)] [Row(0.000000, 3, 3, Double.NegativeInfinity)]
[Row(0.000000, 10, 0, 0.0)] [Row(0.000000, 10, 0, 0.0)]
[Row(0.000000, 10, 1, -inf)] [Row(0.000000, 10, 1, Double.NegativeInfinity)]
[Row(0.000000, 10, 10, -inf)] [Row(0.000000, 10, 10, Double.NegativeInfinity)]
[Row(0.300000, 1, 0, -0.3566749439387324423539544041072745145718090708995)] [Row(0.300000, 1, 0, -0.3566749439387324423539544041072745145718090708995)]
[Row(0.300000, 1, 1, -1.2039728043259360296301803719337238685164245381839)] [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, 0, -1.0700248318161973270618632123218235437154272126985)]
[Row(0.300000, 3, 1, -0.81871040353529122294284394322574719301255212216016)] [Row(0.300000, 3, 1, -0.81871040353529122294284394322574719301255212216016)]
[Row(0.300000, 3, 3, -3.6119184129778080888905411158011716055492736145517)] [Row(0.300000, 3, 3, -3.6119184129778080888905411158011716055492736145517)]
[Row(0.300000, 10, 0, -3.566749439387324423539544041072745145718090708995)] [Row(0.300000, 10, 0, -3.566749439387324423539544041072745145718090708995)]
[Row(0.300000, 10, 1, -2.1114622067804823267977785542148302920616046876506)] [Row(0.300000, 10, 1, -2.1114622067804823267977785542148302920616046876506)]
[Row(0.300000, 10, 10, -12.039728043259360296301803719337238685164245381839)] [Row(0.300000, 10, 10, -12.039728043259360296301803719337238685164245381839)]
[Row(1.000000, 1, 0, -inf)] [Row(1.000000, 1, 0, Double.NegativeInfinity)]
[Row(1.000000, 1, 1, 0.0)]
[Row(1.000000, 1, 1, 0.0)] [Row(1.000000, 1, 1, 0.0)]
[Row(1.000000, 3, 0, -inf)] [Row(1.000000, 3, 0, Double.NegativeInfinity)]
[Row(1.000000, 3, 1, -inf)] [Row(1.000000, 3, 1, Double.NegativeInfinity)]
[Row(1.000000, 3, 3, 0.0)] [Row(1.000000, 3, 3, 0.0)]
[Row(1.000000, 10, 0, -inf)] [Row(1.000000, 10, 0, Double.NegativeInfinity)]
[Row(1.000000, 10, 1, -inf)] [Row(1.000000, 10, 1, Double.NegativeInfinity)]
[Row(1.000000, 10, 10, 0.0)] [Row(1.000000, 10, 10, 0.0)]
public void ValidateProbabilityLn(double p, int n, int x, double dln) public void ValidateProbabilityLn(double p, int n, int x, double dln)
{ {
var b = new Binomial(p,n); var b = new Binomial(p,n);
AssertEx.AreEqual(dln, b.ProbabilityLn(x)); AssertHelpers.AlmostEqual(dln, b.ProbabilityLn(x), 14);
} }
[Test] [Test]

16
src/UnitTests/DistributionTests/Discrete/CategoricalTests.cs

@ -54,7 +54,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
[Test] [Test]
public void CanCreateCategorical() public void CanCreateCategorical()
{ {
var m = new Categorical(largeP, 4); var m = new Categorical(largeP);
AssertEx.AreEqual<double[]>(largeP, m.P); AssertEx.AreEqual<double[]>(largeP, m.P);
} }
@ -62,14 +62,14 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
[ExpectedException(typeof(ArgumentOutOfRangeException))] [ExpectedException(typeof(ArgumentOutOfRangeException))]
public void CategoricalCreateFailsWithNegativeRatios() public void CategoricalCreateFailsWithNegativeRatios()
{ {
var m = new Categorical(badP, 4); var m = new Categorical(badP);
} }
[Test] [Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))] [ExpectedException(typeof(ArgumentOutOfRangeException))]
public void CategoricalCreateFailsWithAllZeroRatios() public void CategoricalCreateFailsWithAllZeroRatios()
{ {
var m = new Categorical(badP2, 4); var m = new Categorical(badP2);
} }
[Test] [Test]
@ -82,7 +82,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
[Test] [Test]
public void CanSetProbability() public void CanSetProbability()
{ {
var b = new Categorical(largeP, 4); var b = new Categorical(largeP);
b.P = smallP; b.P = smallP;
} }
@ -90,27 +90,27 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
[ExpectedException(typeof(ArgumentOutOfRangeException))] [ExpectedException(typeof(ArgumentOutOfRangeException))]
public void SetProbabilityFails() public void SetProbabilityFails()
{ {
var b = new Categorical(largeP, 4); var b = new Categorical(largeP);
b.P = badP; b.P = badP;
} }
[Test] [Test]
public void CanSampleStatic() public void CanSampleStatic()
{ {
var d = Categorical.Sample(new Random(), largeP, 4); var d = Categorical.Sample(new Random(), largeP);
} }
[Test] [Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))] [ExpectedException(typeof(ArgumentOutOfRangeException))]
public void FailSampleStatic() public void FailSampleStatic()
{ {
var d = Categorical.Sample(new Random(), badP, 4); var d = Categorical.Sample(new Random(), badP);
} }
[Test] [Test]
public void CanSample() public void CanSample()
{ {
var n = new Categorical(largeP, 4); var n = new Categorical(largeP);
var d = n.Sample(); var d = n.Sample();
} }
} }

2
src/UnitTests/DistributionTests/Multivariate/MultinomialTests.cs

@ -76,7 +76,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
public void ValidateToString() public void ValidateToString()
{ {
var b = new Multinomial(smallP, 4); 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] [Test]

Loading…
Cancel
Save