diff --git a/src/Numerics/Distributions/Categorical.cs b/src/Numerics/Distributions/Categorical.cs index 538b72af..4394eba2 100644 --- a/src/Numerics/Distributions/Categorical.cs +++ b/src/Numerics/Distributions/Categorical.cs @@ -48,6 +48,9 @@ namespace MathNet.Numerics.Distributions /// does not have to be normalized and sum to 1. The reason is that some vectors can't be exactly normalized /// to sum to 1 in floating point representation. /// + /// + /// Support: 0..k where k = length(probability mass array)-1 + /// public class Categorical : IDiscreteDistribution { System.Random _random; @@ -215,14 +218,14 @@ namespace MathNet.Numerics.Distributions { get { + // Mean = E[X] = Sum(x * p(x), x=0..N-1) + // where f(x) is the probability mass function, and N is the number of categories. + var sum = 0.0; - // Mean = Sum(x * f(x), x=0..N) - // where f(x) is the probability mass function, and N is the maximum value. for (int i = 0; i < _pmfNormalized.Length; i++) { sum += i * _pmfNormalized[i]; } - return sum; } } @@ -232,7 +235,7 @@ namespace MathNet.Numerics.Distributions /// public double StdDev { - get { return _pmfNormalized.StandardDeviation(); } + get { return Math.Sqrt(Variance); } } /// @@ -240,7 +243,18 @@ namespace MathNet.Numerics.Distributions /// public double Variance { - get { return _pmfNormalized.Variance(); } + get + { + // Variance = E[(X-E[X])^2] = E[X^2] - (E[X])^2 = Sum(p(x) * (x - E[X])^2), x=0..N-1) + var m = Mean; + var sum = 0.0; + for (int i = 0; i < _pmfNormalized.Length; i++) + { + var r = i - m; + sum += r*r*_pmfNormalized[i]; + } + return sum; + } } /// @@ -290,7 +304,7 @@ namespace MathNet.Numerics.Distributions /// public double Median { - get { return _pmfNormalized.Median(); } + get { return InverseCumulativeDistribution(0.5); } } /// diff --git a/src/UnitTests/DistributionTests/Discrete/CategoricalTests.cs b/src/UnitTests/DistributionTests/Discrete/CategoricalTests.cs index 70fc55ca..9b805c24 100644 --- a/src/UnitTests/DistributionTests/Discrete/CategoricalTests.cs +++ b/src/UnitTests/DistributionTests/Discrete/CategoricalTests.cs @@ -171,8 +171,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Discrete [TestCase(new double[] { 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 }, 5)] public void ValidateMean(double[] p, double mean) { - var n = new Categorical(p); - AssertHelpers.AlmostEqual(mean, n.Mean, 14); + Assert.That(new Categorical(p).Mean, Is.EqualTo(mean).Within(1e-14)); } /// @@ -187,8 +186,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Discrete [TestCase(new double[] { 1, 0, 1 }, 1)] public void ValidateStdDev(double[] p, double stdDev) { - var n = new Categorical(p); - AssertHelpers.AlmostEqual(stdDev, n.StdDev, 14); + Assert.That(new Categorical(p).StdDev, Is.EqualTo(stdDev).Within(1e-14)); } /// @@ -203,8 +201,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Discrete [TestCase(new double[] { 1, 0, 1 }, 1)] public void ValidateVariance(double[] p, double variance) { - var n = new Categorical(p); - AssertHelpers.AlmostEqual(variance, n.Variance, 14); + Assert.That(new Categorical(p).Variance, Is.EqualTo(variance).Within(1e-14)); } /// @@ -219,13 +216,12 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Discrete // P(X < 5) = (1+2+6+3+2)/29 = 14/29 < 0.5. // P(X <= 5) = 19/29 > 0.5. [TestCase(new double[] { 1, 2, 6, 3, 2, 5, 1, 1, 0, 1, 7 }, 5)] - // TODO: Find out the expected behavour of Median in ambiguous cases like the following: + // TODO: Find out the expected behavior of Median in ambiguous cases like the following: //[TestCase(new double[] { 0, 0.5, 0.5 }, ???)] //[TestCase(new double[] { 1, 0, 1 }, ???)] public void ValidateMedian(double[] p, int median) { - var n = new Categorical(p); - Assert.AreEqual(median, n.Median); + Assert.That(new Categorical(p).Median, Is.EqualTo(median).Within(1e-14)); } ///