Browse Source

Added Gamma & Beta distribution unit tests.

Signed-off-by: Christoph Ruegg <git@cdrnet.ch>
pull/2/head
jvangael 17 years ago
committed by Christoph Ruegg
parent
commit
453023c35b
  1. 83
      src/Managed.UnitTests/DistributionTests/CommonDistributionTests.cs
  2. 378
      src/Managed.UnitTests/DistributionTests/Continuous/BetaTests.cs
  3. 21
      src/Managed.UnitTests/DistributionTests/Continuous/ContinuousUniformTests.cs
  4. 393
      src/Managed.UnitTests/DistributionTests/Continuous/GammaTests.cs
  5. 21
      src/Managed.UnitTests/DistributionTests/Continuous/NormalTests.cs
  6. 3
      src/Managed.UnitTests/Managed.UnitTests.csproj
  7. 324
      src/Managed/Distributions/Continuous/Beta.cs
  8. 504
      src/Managed/Distributions/Continuous/Gamma.cs
  9. 2
      src/Managed/Managed.csproj
  10. 32
      src/Managed/SpecialFunctions.cs
  11. 6
      src/Native.UnitTests/Native.UnitTests.csproj
  12. 6
      src/Native/Native.csproj

83
src/Managed.UnitTests/DistributionTests/CommonDistributionTests.cs

@ -0,0 +1,83 @@
// <copyright file="GammaTests.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://mathnet.opensourcedotnet.info
//
// Copyright (c) 2009 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.UnitTests.DistributionTests
{
using System;
using System.Linq;
using MbUnit.Framework;
using MathNet.Numerics.Distributions;
[TestFixture]
public class CommonDistributionTests
{
private IDistribution[] dists;
[SetUp]
public void SetupDistributions()
{
dists = new IDistribution[4];
dists[0] = new Beta(1.0, 1.0);
dists[1] = new ContinuousUniform(0.0, 1.0);
dists[2] = new Gamma(1.0, 1.0);
dists[3] = new Normal(0.0, 1.0);
}
[Test]
[Row(0)]
[Row(1)]
[Row(2)]
[Row(3)]
public void CanCreateNormal(int i)
{
Assert.IsNotNull(dists[i].RandomSource);
}
[Test]
[Row(0)]
[Row(1)]
[Row(2)]
[Row(3)]
public void CanSetRandomSource(int i)
{
dists[i].RandomSource = new Random();
}
[Test]
[Row(0)]
[Row(1)]
[Row(2)]
[Row(3)]
[ExpectedException(typeof(ArgumentNullException))]
public void FailSetRandomSourceWithNullReference(int i)
{
dists[i].RandomSource = null;
}
}
}

378
src/Managed.UnitTests/DistributionTests/Continuous/BetaTests.cs

@ -0,0 +1,378 @@
// <copyright file="BetaTests.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://mathnet.opensourcedotnet.info
//
// Copyright (c) 2009 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.UnitTests.DistributionTests
{
using System;
using System.Linq;
using MbUnit.Framework;
using MathNet.Numerics.Distributions;
[TestFixture]
public class BetaTests
{
[SetUp]
public void SetUp()
{
Control.CheckDistributionParameters = true;
}
[Test, MultipleAsserts]
[Row(0.0, 0.0)]
[Row(0.0, 0.1)]
[Row(1.0, 0.0)]
[Row(1.0, 1.0)]
[Row(9.0, 1.0)]
[Row(5.0, 100.0)]
[Row(1.0, Double.PositiveInfinity)]
[Row(Double.PositiveInfinity, 1.0)]
[Row(0.0, Double.PositiveInfinity)]
[Row(Double.PositiveInfinity, 0.0)]
public void CanCreateBeta(double a, double b)
{
var n = new Beta(a, b);
AssertEx.AreEqual<double>(a, n.A);
AssertEx.AreEqual<double>(b, n.B);
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
[Row(Double.NaN, 1.0)]
[Row(1.0, Double.NaN)]
[Row(Double.NaN, Double.NaN)]
[Row(1.0, -1.0)]
[Row(-1.0, 1.0)]
[Row(-1.0, -1.0)]
public void BetaCreateFailsWithBadParameters(double a, double b)
{
var n = new Beta(a, b);
}
[Test]
public void ValidateToString()
{
var n = new Beta(1.0, 2.0);
AssertEx.AreEqual<string>("Beta(A = 1, B = 2)", n.ToString());
}
[Test]
[Row(-0.0)]
[Row(0.0)]
[Row(0.1)]
[Row(1.0)]
[Row(10.0)]
[Row(Double.PositiveInfinity)]
public void CanSetShapeA(double a)
{
var n = new Beta(1.0, 1.0);
n.A = a;
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void SetShapeAFailsWithNegativeA()
{
var n = new Beta(1.0, 1.0);
n.A = -1.0;
}
[Test]
[Row(-0.0)]
[Row(0.0)]
[Row(0.1)]
[Row(1.0)]
[Row(10.0)]
[Row(Double.PositiveInfinity)]
public void CanSetShapeB(double b)
{
var n = new Beta(1.0, 1.0);
n.B = b;
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void SetShapeBFailsWithNegativeB()
{
var n = new Beta(1.0, 1.0);
n.B = -1.0;
}
[Test]
[Row(0.0, 0.0, 0.5)]
[Row(0.0, 0.1, 0.1)]
[Row(1.0, 0.0, 1.0)]
[Row(1.0, 1.0, 0.5)]
[Row(9.0, 1.0, 0.9)]
[Row(5.0, 100.0, 0.047619047619047619047616)]
[Row(1.0, Double.PositiveInfinity, 1.0)]
[Row(Double.PositiveInfinity, 1.0, 0.0)]
[Row(0.0, Double.PositiveInfinity, 1.0)]
[Row(Double.PositiveInfinity, 0.0, 0.0)]
public void ValidateMean(double a, double b, double mean)
{
var n = new Beta(a, b);
AssertEx.AreEqual<double>(mean, n.Mean);
}
[Test]
[Row(0.0, 0.0, 0.5)]
[Row(0.0, 0.1, 0.1)]
[Row(1.0, 0.0, 1.0)]
[Row(1.0, 1.0, 0.0)]
[Row(9.0, 1.0, -1.3083356884473304939016015849561625204060922267565917)]
[Row(5.0, 100.0, -2.5201623187602743679459255108827601222133603091753153)]
[Row(1.0, Double.PositiveInfinity, 0.0)]
[Row(Double.PositiveInfinity, 1.0, 0.0)]
[Row(0.0, Double.PositiveInfinity, 0.0)]
[Row(Double.PositiveInfinity, 0.0, 0.0)]
public void ValidateEntropy(double a, double b, double entropy)
{
var n = new Beta(a, b);
AssertEx.AreEqual<double>(entropy, n.Entropy);
}
[Test]
[Row(0.0, 0.0, 0.0)]
[Row(0.0, 0.1, 2.0)]
[Row(1.0, 0.0, -2.0)]
[Row(1.0, 1.0, 0.0)]
[Row(9.0, 1.0, -1.4740554623801777107177478829647496373009282424841579)]
[Row(5.0, 100.0, 0.81759410927553430354583159143895018978562196953345572)]
[Row(1.0, Double.PositiveInfinity, 2.0)]
[Row(Double.PositiveInfinity, 1.0, -2.0)]
[Row(0.0, Double.PositiveInfinity, 2.0)]
[Row(Double.PositiveInfinity, 0.0, -2.0)]
public void ValidateSkewness(double a, double b, double skewness)
{
var n = new Beta(a, b);
AssertEx.AreEqual<double>(skewness, n.Skewness);
}
[Test]
[Row(0.0, 0.0, 0.5)]
[Row(0.0, 0.1, 1.0)]
[Row(1.0, 0.0, 0.0)]
[Row(1.0, 1.0, 0.5)]
[Row(9.0, 1.0, 1.0)]
[Row(5.0, 100.0, 0.038834951456310676243255386452801758423447608947753906)]
[Row(1.0, Double.PositiveInfinity, 0.0)]
[Row(Double.PositiveInfinity, 1.0, 1.0)]
[Row(0.0, Double.PositiveInfinity, 0.0)]
[Row(Double.PositiveInfinity, 0.0, 1.0)]
public void ValidateMode(double a, double b, double mode)
{
var n = new Beta(a, b);
AssertEx.AreEqual<double>(mode, n.Mode);
}
[Test]
[ExpectedException(typeof(NotSupportedException))]
[Row(0.0, 0.0)]
[Row(0.0, 0.1)]
[Row(1.0, 0.0)]
[Row(1.0, 1.0)]
[Row(9.0, 1.0)]
[Row(5.0, 100.0)]
[Row(1.0, Double.PositiveInfinity)]
[Row(Double.PositiveInfinity, 1.0)]
[Row(0.0, Double.PositiveInfinity)]
[Row(Double.PositiveInfinity, 0.0)]
public void ValidateMedian(double a, double b)
{
var n = new Beta(a, 1.0);
var m = n.Median;
}
[Test]
public void ValidateMinimum()
{
var n = new Beta(1.0, 1.0);
AssertEx.AreEqual<double>(0.0, n.Minimum);
}
[Test]
public void ValidateMaximum()
{
var n = new Beta(1.0, 1.0);
AssertEx.AreEqual<double>(1.0, n.Maximum);
}
[Test]
public void CanSampleStatic()
{
var d = Beta.Sample(new Random(), 2.0, 3.0);
}
[Test]
public void CanSampleSequenceStatic()
{
var ied = Beta.Samples(new Random(), 2.0, 3.0);
var arr = ied.Take(5).ToArray();
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void FailSampleStatic()
{
var d = Beta.Sample(new Random(), 1.0, -1.0);
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void FailSampleSequenceStatic()
{
var ied = Beta.Samples(new Random(), 1.0, -1.0).First();
}
[Test]
public void CanSample()
{
var n = new Beta(2.0, 3.0);
var d = n.Sample();
}
[Test]
public void CanSampleSequence()
{
var n = new Beta(2.0, 3.0);
var ied = n.Samples();
var e = ied.Take(5).ToArray();
}
[Test]
[Row(0.0, 0.0, 0.0, Double.PositiveInfinity)]
[Row(0.0, 0.0, 0.5, 0.0)]
[Row(0.0, 0.0, 1.0, Double.PositiveInfinity)]
[Row(0.0, 0.1, 0.0, Double.PositiveInfinity)]
[Row(0.0, 0.1, 0.5, 0.0)]
[Row(0.0, 0.1, 1.0, 0.0)]
[Row(1.0, 0.0, 0.0, 0.0)]
[Row(1.0, 0.0, 0.5, 0.0)]
[Row(1.0, 0.0, 1.0, Double.PositiveInfinity)]
[Row(1.0, 1.0, 0.0, 1.0)]
[Row(1.0, 1.0, 0.5, 1.0)]
[Row(1.0, 1.0, 1.0, 1.0)]
[Row(9.0, 1.0, 0.0, 0.0)]
[Row(9.0, 1.0, 0.5, 0.035155378090821160189479427593561667617370600927556366)]
[Row(9.0, 1.0, 1.0, 8.9997767912502170085067334639517869100468738374544298)]
[Row(5.0, 100.0, 0.0, 0.0)]
[Row(5.0, 100.0, 0.5, 1.0881845516040810386311829462908430145307026037926335e-21)]
[Row(5.0, 100.0, 1.0, 0.0)]
[Row(1.0, Double.PositiveInfinity, 0.0, Double.PositiveInfinity)]
[Row(1.0, Double.PositiveInfinity, 0.5, 0.0)]
[Row(1.0, Double.PositiveInfinity, 1.0, 0.0)]
[Row(Double.PositiveInfinity, 1.0, 0.0, 0.0)]
[Row(Double.PositiveInfinity, 1.0, 0.5, 0.0)]
[Row(Double.PositiveInfinity, 1.0, 1.0, Double.PositiveInfinity)]
[Row(0.0, Double.PositiveInfinity, 0.0, Double.PositiveInfinity)]
[Row(0.0, Double.PositiveInfinity, 0.5, 0.0)]
[Row(0.0, Double.PositiveInfinity, 1.0, 0.0)]
[Row(Double.PositiveInfinity, 0.0, 0.0, 0.0)]
[Row(Double.PositiveInfinity, 0.0, 0.5, 0.0)]
[Row(Double.PositiveInfinity, 0.0, 1.0, Double.PositiveInfinity)]
public void ValidateDensity(double a, double b, double x, double pdf)
{
var n = new Beta(a, b);
AssertHelpers.AlmostEqual(pdf, n.Density(x), 15);
}
[Test]
[Row(0.0, 0.0, 0.0, Double.PositiveInfinity)]
[Row(0.0, 0.0, 0.5, Double.NegativeInfinity)]
[Row(0.0, 0.0, 1.0, Double.PositiveInfinity)]
[Row(0.0, 0.1, 0.0, Double.PositiveInfinity)]
[Row(0.0, 0.1, 0.5, Double.NegativeInfinity)]
[Row(0.0, 0.1, 1.0, Double.NegativeInfinity)]
[Row(1.0, 0.0, 0.0, Double.NegativeInfinity)]
[Row(1.0, 0.0, 0.5, Double.NegativeInfinity)]
[Row(1.0, 0.0, 1.0, Double.PositiveInfinity)]
[Row(1.0, 1.0, 0.0, 0.0)]
[Row(1.0, 1.0, 0.5, 0.0)]
[Row(1.0, 1.0, 1.0, 0.0)]
[Row(9.0, 1.0, 0.0, Double.NegativeInfinity)]
[Row(9.0, 1.0, 0.5, -3.3479528671433430925473664978203611353090199592365404)]
[Row(9.0, 1.0, 1.0, 2.197199776056471990627201569939718235188380979206923)]
[Row(5.0, 100.0, 0.0, Double.NegativeInfinity)]
[Row(5.0, 100.0, 0.5, -51.447830024537682154565870837960406410586196074573801)]
[Row(5.0, 100.0, 1.0, Double.NegativeInfinity)]
[Row(1.0, Double.PositiveInfinity, 0.0, Double.PositiveInfinity)]
[Row(1.0, Double.PositiveInfinity, 0.5, Double.NegativeInfinity)]
[Row(1.0, Double.PositiveInfinity, 1.0, Double.NegativeInfinity)]
[Row(Double.PositiveInfinity, 1.0, 0.0, Double.NegativeInfinity)]
[Row(Double.PositiveInfinity, 1.0, 0.5, Double.NegativeInfinity)]
[Row(Double.PositiveInfinity, 1.0, 1.0, Double.PositiveInfinity)]
[Row(0.0, Double.PositiveInfinity, 0.0, Double.PositiveInfinity)]
[Row(0.0, Double.PositiveInfinity, 0.5, Double.NegativeInfinity)]
[Row(0.0, Double.PositiveInfinity, 1.0, Double.NegativeInfinity)]
[Row(Double.PositiveInfinity, 0.0, 0.0, Double.NegativeInfinity)]
[Row(Double.PositiveInfinity, 0.0, 0.5, Double.NegativeInfinity)]
[Row(Double.PositiveInfinity, 0.0, 1.0, Double.PositiveInfinity)]
public void ValidateDensityLn(double a, double b, double x, double pdfln)
{
var n = new Beta(a, b);
AssertHelpers.AlmostEqual(pdfln, n.DensityLn(x), 15);
}
[Test]
[Row(0.0, 0.0, 0.0, 0.5)]
[Row(0.0, 0.0, 0.5, 0.5)]
[Row(0.0, 0.0, 1.0, 1.0)]
[Row(0.0, 0.1, 0.0, 1.0)]
[Row(0.0, 0.1, 0.5, 1.0)]
[Row(0.0, 0.1, 1.0, 1.0)]
[Row(1.0, 0.0, 0.0, 0.0)]
[Row(1.0, 0.0, 0.5, 0.0)]
[Row(1.0, 0.0, 1.0, 1.0)]
[Row(1.0, 1.0, 0.0, 0.0)]
[Row(1.0, 1.0, 0.5, 0.5)]
[Row(1.0, 1.0, 1.0, 1.0)]
[Row(9.0, 1.0, 0.0, 0.0)]
[Row(9.0, 1.0, 0.5, 0.00195313)]
[Row(9.0, 1.0, 1.0, 1.0)]
[Row(5.0, 100.0, 0.0, 0.0)]
[Row(5.0, 100.0, 0.5, 1.0)]
[Row(5.0, 100.0, 1.0, 1.0)]
[Row(1.0, Double.PositiveInfinity, 0.0, 1.0)]
[Row(1.0, Double.PositiveInfinity, 0.5, 1.0)]
[Row(1.0, Double.PositiveInfinity, 1.0, 1.0)]
[Row(Double.PositiveInfinity, 1.0, 0.0, 0.0)]
[Row(Double.PositiveInfinity, 1.0, 0.5, 0.0)]
[Row(Double.PositiveInfinity, 1.0, 1.0, 1.0)]
[Row(0.0, Double.PositiveInfinity, 0.0, 1.0)]
[Row(0.0, Double.PositiveInfinity, 0.5, 1.0)]
[Row(0.0, Double.PositiveInfinity, 1.0, 1.0)]
[Row(Double.PositiveInfinity, 0.0, 0.0, 0.0)]
[Row(Double.PositiveInfinity, 0.0, 0.5, 0.0)]
[Row(Double.PositiveInfinity, 0.0, 1.0, 1.0)]
public void ValidateCumulativeDistribution(double a, double b, double x, double cdf)
{
var n = new Beta(a, b);
AssertHelpers.AlmostEqual(cdf, n.CumulativeDistribution(x), 15);
}
}
}

21
src/Managed.UnitTests/DistributionTests/Continuous/ContinuousUniformTests.cs

@ -83,21 +83,6 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
AssertEx.AreEqual<string>("ContinuousUniform(Lower = 1, Upper = 2)", n.ToString());
}
[Test]
public void CanGetRandomSource()
{
var n = new ContinuousUniform();
var rs = n.RandomSource;
Assert.IsNotNull(rs);
}
[Test]
public void CanSetRandomSource()
{
var n = new ContinuousUniform();
n.RandomSource = new Random();
}
[Test]
[Row(-10.0)]
[Row(-0.0)]
@ -272,11 +257,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
public void CanSampleSequenceStatic()
{
var ied = ContinuousUniform.Samples(new Random(), 0.0, 1.0);
var e = ied.GetEnumerator();
e.MoveNext();
var d = e.Current;
e.MoveNext();
var g = e.Current;
var arr = ied.Take(5).ToArray();
}
[Test]

393
src/Managed.UnitTests/DistributionTests/Continuous/GammaTests.cs

@ -0,0 +1,393 @@
// <copyright file="GammaTests.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://mathnet.opensourcedotnet.info
//
// Copyright (c) 2009 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.UnitTests.DistributionTests
{
using System;
using System.Linq;
using MbUnit.Framework;
using MathNet.Numerics.Distributions;
[TestFixture]
public class GammaTests
{
[SetUp]
public void SetUp()
{
Control.CheckDistributionParameters = true;
}
[Test, MultipleAsserts]
[Row(0.0, 0.0)]
[Row(1.0, 0.1)]
[Row(1.0, 1.0)]
[Row(10.0, 10.0)]
[Row(10.0, 1.0)]
[Row(10.0, Double.PositiveInfinity)]
public void CanCreateGamma(double shape, double invScale)
{
var n = new Gamma(shape, invScale);
AssertEx.AreEqual<double>(shape, n.Shape);
AssertEx.AreEqual<double>(invScale, n.InvScale);
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
[Row(Double.NaN, 1.0)]
[Row(1.0, Double.NaN)]
[Row(Double.NaN, Double.NaN)]
[Row(1.0, -1.0)]
[Row(-1.0, 1.0)]
[Row(-1.0, -1.0)]
public void GammaCreateFailsWithBadParameters(double shape, double invScale)
{
var n = new Gamma(shape, invScale);
}
[Test, MultipleAsserts]
[Row(0.0, 0.0)]
[Row(1.0, 0.1)]
[Row(1.0, 1.0)]
[Row(10.0, 10.0)]
[Row(10.0, 1.0)]
[Row(10.0, Double.PositiveInfinity)]
public void CanCreateGammaWithShapeInvScale(double shape, double invScale)
{
var n = Gamma.WithShapeInvScale(shape, invScale);
AssertEx.AreEqual<double>(shape, n.Shape);
AssertEx.AreEqual<double>(invScale, n.InvScale);
}
[Test, MultipleAsserts]
[Row(0.0, 0.0)]
[Row(1.0, 0.1)]
[Row(1.0, 1.0)]
[Row(10.0, 10.0)]
[Row(10.0, 1.0)]
[Row(10.0, Double.PositiveInfinity)]
public void CanCreateGammaWithShapeScale(double shape, double scale)
{
var n = Gamma.WithShapeScale(shape, scale);
AssertEx.AreEqual<double>(shape, n.Shape);
AssertEx.AreEqual<double>(scale, n.Scale);
}
[Test]
public void ValidateToString()
{
var n = new Gamma(1.0, 2.0);
AssertEx.AreEqual<string>("Gamma(Shape = 1, Inverse Scale = 2)", n.ToString());
}
[Test]
[Row(-0.0)]
[Row(0.0)]
[Row(0.1)]
[Row(1.0)]
[Row(10.0)]
[Row(Double.PositiveInfinity)]
public void CanSetShape(double shape)
{
var n = new Gamma(1.0, 1.0);
n.Shape = shape;
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void SetShapeFailsWithNegativeShape()
{
var n = new Gamma(1.0, 1.0);
n.Shape = -1.0;
}
[Test]
[Row(-0.0)]
[Row(0.0)]
[Row(0.1)]
[Row(1.0)]
[Row(10.0)]
[Row(Double.PositiveInfinity)]
public void CanSetScale(double scale)
{
var n = new Gamma(1.0, 1.0);
n.Scale = scale;
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void SetScaleFailsWithNegativeScale()
{
var n = new Gamma(1.0, 1.0);
n.Scale = -1.0;
}
[Test]
[Row(-0.0)]
[Row(0.0)]
[Row(0.1)]
[Row(1.0)]
[Row(10.0)]
[Row(Double.PositiveInfinity)]
public void CanSetInvScale(double invScale)
{
var n = new Gamma(1.0, 1.0);
n.InvScale = invScale;
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void SetInvScaleFailsWithNegativeInvScale()
{
var n = new Gamma(1.0, 1.0);
n.InvScale = -1.0;
}
[Test]
[Row(0.0, 0.0, 0.0)]
[Row(1.0, 0.1, 10.0)]
[Row(1.0, 1.0, 1.0)]
[Row(10.0, 10.0, 1.0)]
[Row(10.0, 1.0, 10.0)]
[Row(10.0, Double.PositiveInfinity, 0.0)]
public void CanGetMean(double shape, double invScale, double mean)
{
var n = new Gamma(shape, invScale);
AssertEx.AreEqual<double>(mean, n.Mean);
}
[Test]
[Row(0.0, 0.0, 0.0)]
[Row(1.0, 0.1, 100.0)]
[Row(1.0, 1.0, 1.0)]
[Row(10.0, 10.0, 0.1)]
[Row(10.0, 1.0, 10.0)]
[Row(10.0, Double.PositiveInfinity, 0.0)]
public void CanGetVariance(double shape, double invScale, double var)
{
var n = new Gamma(shape, invScale);
AssertEx.AreEqual<double>(var, n.Variance);
}
[Test]
[Row(0.0, 0.0, 0.0)]
[Row(1.0, 0.1, 10.0)]
[Row(1.0, 1.0, 1.0)]
[Row(10.0, 10.0, 0.31622776601683794197697302588502426416723164097476643)]
[Row(10.0, 1.0, 3.1622776601683793319988935444327185337195551393252168)]
[Row(10.0, Double.PositiveInfinity, 0.0)]
public void CanGetStdDev(double shape, double invScale, double sdev)
{
var n = new Gamma(shape, invScale);
AssertHelpers.AlmostEqual(sdev, n.StdDev, 15);
}
[Test]
[Row(0.0, 0.0, Double.PositiveInfinity)]
[Row(1.0, 0.1, 3.3025850929940456285068402234265387271634735938763824)]
[Row(1.0, 1.0, 1.0)]
[Row(10.0, 10.0, 0.23346908548693395836262094490967812177376750477943892)]
[Row(10.0, 1.0, 2.5360541784809796423806123995940423293748689934081866)]
[Row(10.0, Double.PositiveInfinity, 0.0)]
public void ValidateEntropy(double shape, double invScale, double entropy)
{
var n = new Gamma(shape, invScale);
AssertHelpers.AlmostEqual(entropy, n.Entropy, 15);
}
[Test]
[Row(0.0, 0.0, Double.PositiveInfinity)]
[Row(1.0, 0.1, 2.0)]
[Row(1.0, 1.0, 2.0)]
[Row(10.0, 10.0, 0.63245553203367586639977870888654370674391102786504337)]
[Row(10.0, 1.0, 0.63245553203367586639977870888654370674391102786504337)]
[Row(10.0, Double.PositiveInfinity, 0.0)]
public void ValidateSkewness(double shape, double invScale, double skewness)
{
var n = new Gamma(shape, invScale);
AssertHelpers.AlmostEqual(skewness, n.Skewness, 15);
}
[Test]
[Row(0.0, 0.0, Double.PositiveInfinity)]
[Row(1.0, 0.1, 0.0)]
[Row(1.0, 1.0, 0.0)]
[Row(10.0, 10.0, 0.9)]
[Row(10.0, 1.0, 9.0)]
[Row(10.0, Double.PositiveInfinity, 10.0)]
public void ValidateMode(double shape, double invScale, double mode)
{
var n = new Gamma(shape, invScale);
AssertEx.AreEqual<double>(mode, n.Mode);
}
[Test]
[ExpectedException(typeof(NotSupportedException))]
[Row(0.0, 0.0)]
[Row(1.0, 0.1)]
[Row(1.0, 1.0)]
[Row(10.0, 10.0)]
[Row(10.0, 1.0)]
[Row(10.0, Double.PositiveInfinity)]
public void ValidateMedian(double shape, double invScale, double mode)
{
var n = new Gamma(shape, invScale);
var median = n.Median;
}
[Test]
public void ValidateMinimum()
{
var n = new Gamma(1.0,1.0);
AssertEx.AreEqual<double>(0.0, n.Minimum);
}
[Test]
public void ValidateMaximum()
{
var n = new Gamma(1.0, 1.0);
AssertEx.AreEqual<double>(System.Double.PositiveInfinity, n.Maximum);
}
[Test]
[Row(0.0, 0.0, 0.0, 0.0)]
[Row(0.0, 0.0, 1.0, 0.0)]
[Row(0.0, 0.0, 10.0, 0.0)]
[Row(1.0, 0.1, 0.0, 0.10000000000000000555111512312578270211815834045410156)]
[Row(1.0, 0.1, 1.0, 0.099004983374916810660915381324116279472953858297127391)]
[Row(1.0, 0.1, 10.0, 0.036787944117144234201693506390001264039984687455876246)]
[Row(1.0, 1.0, 0.0, 1.0)]
[Row(1.0, 1.0, 1.0, 0.36787944117144232159552377016146086744581113103176804)]
[Row(1.0, 1.0, 10.0, 0.000045399929762484851535591515560550610237918088866564953)]
[Row(10.0, 10.0, 0.0, 0.0)]
[Row(10.0, 10.0, 1.0, 1.2511003572113329898476497894772544708420990097708588)]
[Row(10.0, 10.0, 10.0, 1.0251532120868705806216092933926141802686541811003037e-30)]
[Row(10.0, 1.0, 0.0, 0.0)]
[Row(10.0, 1.0, 1.0, 0.0000010137771196302974029859010421116095333052555418644397)]
[Row(10.0, 1.0, 10.0, 0.12511003572113329898476497894772544708420990097708601)]
[Row(10.0, Double.PositiveInfinity, 0.0, 0.0)]
[Row(10.0, Double.PositiveInfinity, 1.0, 0.0)]
[Row(10.0, Double.PositiveInfinity, 10.0, Double.PositiveInfinity)]
public void ValidateDensity(double shape, double invScale, double x, double pdf)
{
var n = new Gamma(shape, invScale);
AssertHelpers.AlmostEqual(pdf, n.Density(x), 15);
}
[Test]
[Row(0.0, 0.0, 0.0, Double.NegativeInfinity)]
[Row(0.0, 0.0, 1.0, Double.NegativeInfinity)]
[Row(0.0, 0.0, 10.0, Double.NegativeInfinity)]
[Row(1.0, 0.1, 0.0, -2.3025850929940456285068402234265387271634735938763824)]
[Row(1.0, 0.1, 1.0, -2.312585092994045630449730516520562672904829013035318)]
[Row(1.0, 0.1, 10.0, -3.3025850929940456285068402234265387271634735938763824)]
[Row(1.0, 1.0, 0.0, 0.0)]
[Row(1.0, 1.0, 1.0, -1.0)]
[Row(1.0, 1.0, 10.0, -10.0)]
[Row(10.0, 10.0, 0.0, Double.NegativeInfinity)]
[Row(10.0, 10.0, 1.0, 0.22402344985898722897219667227693591172986563062456522)]
[Row(10.0, 10.0, 10.0, -69.052710713194601614865880235563786219860220971716511)]
[Row(10.0, 1.0, 0.0, Double.NegativeInfinity)]
[Row(10.0, 1.0, 1.0, -13.801827480081469611207717874566706164281149255663166)]
[Row(10.0, 1.0, 10.0, -2.0785616431350584550457947824074282958712358580042068)]
[Row(10.0, Double.PositiveInfinity, 0.0, Double.NegativeInfinity)]
[Row(10.0, Double.PositiveInfinity, 1.0, Double.NegativeInfinity)]
[Row(10.0, Double.PositiveInfinity, 10.0, Double.PositiveInfinity)]
public void ValidateDensityLn(double shape, double invScale, double x, double pdfln)
{
var n = new Gamma(shape, invScale);
AssertHelpers.AlmostEqual(pdfln, n.DensityLn(x), 15);
}
[Test]
public void CanSampleStatic()
{
var d = Gamma.Sample(new Random(), 1.0, 1.0);
}
[Test]
public void CanSampleSequenceStatic()
{
var ied = Gamma.Samples(new Random(), 1.0, 1.0);
var arr = ied.Take(5).ToArray();
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void FailSampleStatic()
{
var d = Normal.Sample(new Random(), 1.0, -1.0);
}
[Test]
[ExpectedException(typeof(ArgumentOutOfRangeException))]
public void FailSampleSequenceStatic()
{
var ied = Normal.Samples(new Random(), 1.0, -1.0).First();
}
[Test]
public void CanSample()
{
var n = new Normal();
var d = n.Sample();
}
[Test]
public void CanSampleSequence()
{
var n = new Normal();
var ied = n.Samples();
var e = ied.Take(5).ToArray();
}
[Test]
[Row(0.0, 0.0, 0.0, 0.0)]
[Row(0.0, 0.0, 1.0, 0.0)]
[Row(0.0, 0.0, 10.0, 0.0)]
[Row(1.0, 0.1, 0.0, 0.0)]
[Row(1.0, 0.1, 1.0, 0.095162581964040431858607615783064404690935346242622848)]
[Row(1.0, 0.1, 10.0, 0.63212055882855767840447622983853913255418886896823196)]
[Row(1.0, 1.0, 0.0, 0.0)]
[Row(1.0, 1.0, 1.0, 0.63212055882855767840447622983853913255418886896823196)]
[Row(1.0, 1.0, 10.0, 0.99995460007023751514846440848443944938976208191113396)]
[Row(10.0, 10.0, 0.0, 0.0)]
[Row(10.0, 10.0, 1.0, 0.54207028552814779168583514294066541824736464003242184)]
[Row(10.0, 10.0, 10.0, 0.99999999999999999999999999999988746526039157266114706)]
[Row(10.0, 1.0, 0.0, 0.0)]
[Row(10.0, 1.0, 1.0, 0.00000011142547833872067735305068724025236288094949815466035)]
[Row(10.0, 1.0, 10.0, 0.54207028552814779168583514294066541824736464003242184)]
[Row(10.0, Double.PositiveInfinity, 0.0, 0.0)]
[Row(10.0, Double.PositiveInfinity, 1.0, 0.0)]
[Row(10.0, Double.PositiveInfinity, 10.0, 1.0)]
public void ValidateCumulativeDistribution(double shape, double invScale, double x, double cdf)
{
var n = new Gamma(shape, invScale);
AssertHelpers.AlmostEqual(cdf, n.CumulativeDistribution(x), 15);
}
}
}

21
src/Managed.UnitTests/DistributionTests/Continuous/NormalTests.cs

@ -128,21 +128,6 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
AssertEx.AreEqual<string>("Normal(Mean = 1, StdDev = 2)", n.ToString());
}
[Test]
public void CanGetRandomSource()
{
var n = new Normal();
var rs = n.RandomSource;
Assert.IsNotNull(rs);
}
[Test]
public void CanSetRandomSource()
{
var n = new Normal();
n.RandomSource = new Random();
}
[Test]
[Row(-0.0)]
[Row(0.0)]
@ -338,11 +323,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
public void CanSampleSequenceStatic()
{
var ied = Normal.Samples(new Random(), 0.0, 1.0);
var e = ied.GetEnumerator();
e.MoveNext();
var d = e.Current;
e.MoveNext();
var g = e.Current;
var arr = ied.Take(5).ToArray();
}
[Test]

3
src/Managed.UnitTests/Managed.UnitTests.csproj

@ -62,7 +62,10 @@
<Compile Include="AssertHelpers.cs" />
<Compile Include="CombinatoricsTests\CombinatoricsCountingTest.cs" />
<Compile Include="ComplexTests\ComplexTest.cs" />
<Compile Include="DistributionTests\CommonDistributionTests.cs" />
<Compile Include="DistributionTests\Continuous\BetaTests.cs" />
<Compile Include="DistributionTests\Continuous\ContinuousUniformTests.cs" />
<Compile Include="DistributionTests\Continuous\GammaTests.cs" />
<Compile Include="DistributionTests\Continuous\NormalTests.cs" />
<Compile Include="IntegralTransformsTests\DftTest.cs" />
<Compile Include="IntegrationTests\IntegrationTest.cs" />

324
src/Managed/Distributions/Continuous/Beta.cs

@ -0,0 +1,324 @@
// <copyright file="Beta.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://mathnet.opensourcedotnet.info
//
// Copyright (c) 2009 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.Distributions
{
using System;
using System.Collections.Generic;
using Properties;
/// <summary>
/// Implements the Beta distribution. For details about this distribution, see
/// <a href="http://en.wikipedia.org/wiki/Beta_distribution">Wikipedia - Beta distribution</a>.
/// </summary>
/// <remarks><para>The distribution will use the <see cref="System.Random"/> by default.
/// Users can get/set the random number generator by using the <see cref="RandomSource"/> property.</para>
/// <para>The statistics classes will check all the incoming parameters whether they are in the allowed
/// range. This might involve heavy computation. Optionally, by setting Control.CheckDistributionParameters
/// to false, all parameter checks can be turned off.</para></remarks>
public class Beta : IContinuousDistribution
{
/// <summary>
/// Beta shape parameter a.
/// </summary>
private double _shapeA;
/// <summary>
/// Beta shape parameter b.
/// </summary>
private double _shapeB;
/// <summary>
/// Initializes a new instance of the Beta distribution.
/// </summary>
/// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">If any of the Beta parameters are negative.</exception>
public Beta(double a, double b)
{
SetParameters(a, b);
RandomSource = new Random();
}
/// <summary>
/// A string representation of the distribution.
/// </summary>
public override string ToString()
{
return "Beta(A = " + _shapeA + ", B = " + _shapeB + ")";
}
/// <summary>
/// Checks whether the parameters of the distribution are valid.
/// </summary>
/// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <returns>True when the parameters are valid, false otherwise.</returns>
private static bool IsValidParameterSet(double a, double b)
{
if (a < 0.0 || b < 0.0)
{
return false;
}
return true;
}
/// <summary>
/// Sets the parameters of the distribution after checking their validity.
/// </summary>
/// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double a, double b)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
_shapeA = a;
_shapeB = b;
}
/// <summary>
/// Gets or sets the A shape parameter of the Beta distribution.
/// </summary>
public double A
{
get { return _shapeA; }
set { SetParameters(value, _shapeB); }
}
/// <summary>
/// Gets or sets the B shape parameter of the Beta distribution.
/// </summary>
public double B
{
get { return _shapeB; }
set { SetParameters(_shapeA, value); }
}
#region IDistribution implementation
/// <summary>
/// Gets or sets the random number generator which is used to draw random samples.
/// </summary>
public Random RandomSource { get; set; }
/// <summary>
/// Gets the mean of the Beta distribution.
/// </summary>
public double Mean
{
get { return _shapeA / (_shapeA + _shapeB); }
}
/// <summary>
/// Gets the variance of the Beta distribution.
/// </summary>
public double Variance
{
get { return (_shapeA * _shapeB) / ((_shapeA + _shapeB) * (_shapeA + _shapeB) * (_shapeA + _shapeB + 1.0)); }
}
/// <summary>
/// Gets the standard deviation of the Beta distribution.
/// </summary>
public double StdDev
{
get { return Math.Sqrt((_shapeA * _shapeB) / ((_shapeA + _shapeB) * (_shapeA + _shapeB) * (_shapeA + _shapeB + 1.0))); }
}
/// <summary>
/// Gets the entropy of the Beta distribution.
/// </summary>
public double Entropy
{
get
{
return SpecialFunctions.BetaLn(_shapeA, _shapeB)
- (_shapeA - 1.0) * SpecialFunctions.DiGamma(_shapeA)
- (_shapeB - 1.0) * SpecialFunctions.DiGamma(_shapeB)
+ (_shapeA + _shapeB - 2.0) * SpecialFunctions.DiGamma(_shapeA + _shapeB);
}
}
/// <summary>
/// Gets the skewness of the Beta distribution.
/// </summary>
public double Skewness
{
get
{
return 2.0 * (_shapeB - _shapeA) * Math.Sqrt(_shapeA + _shapeB + 1.0)
/ ((_shapeA + _shapeB + 2.0) * Math.Sqrt(_shapeA * _shapeB));
}
}
#endregion
#region IContinuousDistribution implementation
/// <summary>
/// Gets the mode of the Beta distribution.
/// </summary>
public double Mode
{
get { return (_shapeA - 1) / (_shapeA + _shapeB - 2); }
}
/// <summary>
/// Gets the median of the Beta distribution.
/// </summary>
public double Median
{
get { throw new NotSupportedException(); }
}
/// <summary>
/// Gets the minimum of the Beta distribution.
/// </summary>
public double Minimum
{
get { return 0.0; }
}
/// <summary>
/// Gets the maximum of the Beta distribution.
/// </summary>
public double Maximum
{
get { return 1.0; }
}
/// <summary>
/// Computes the density of the Beta distribution.
/// </summary>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
public double Density(double x)
{
double b = SpecialFunctions.Gamma(_shapeA + _shapeB) / (SpecialFunctions.Gamma(_shapeA) * SpecialFunctions.Gamma(_shapeB));
return b * Math.Pow(x, _shapeA - 1.0) * Math.Pow(1.0 - x, _shapeB - 1.0);
}
/// <summary>
/// Computes the log density of the Beta distribution.
/// </summary>
/// <param name="x">The location at which to compute the log density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
public double DensityLn(double x)
{
double b = SpecialFunctions.GammaLn(_shapeA + _shapeB) - SpecialFunctions.GammaLn(_shapeA) - SpecialFunctions.GammaLn(_shapeB);
return b + (_shapeA - 1.0)*Math.Log(x) + (_shapeB - 1.0)*Math.Log(1.0 - x);
}
/// <summary>
/// Computes the cumulative distribution function of the Beta distribution.
/// </summary>
/// <param name="x">The location at which to compute the cumulative density.</param>
/// <returns>the cumulative density at <paramref name="x"/>.</returns>
public double CumulativeDistribution(double x)
{
return SpecialFunctions.BetaRegularized(_shapeA, _shapeB, x);
}
/// <summary>
/// Generates a sample from the Beta distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return SampleBeta(RandomSource, _shapeA, _shapeB);
}
/// <summary>
/// Generates a sequence of samples from the Beta distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
while (true)
{
yield return SampleBeta(RandomSource, _shapeA, _shapeB);
}
}
#endregion
/// <summary>
/// Generates a sample from the normal distribution using the <i>Box-Muller</i> algorithm.
/// </summary>
/// <param name="rng">The random number generator to use.</param>
/// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rng, double a, double b)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleBeta(rng, a, b);
}
/// <summary>
/// Generates a sequence of samples from the normal distribution using the <i>Box-Muller</i> algorithm.
/// </summary>
/// <param name="rng">The random number generator to use.</param>
/// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rng, double a, double b)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleBeta(rng, a, b);
}
}
/// <summary>
/// Samples Beta distributed random variables by sampling two Gamma variables and normalizing.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="a">The A shape parameter.</param>
/// <param name="b">The B shape parameter.</param>
/// <returns>a random number from the Beta distribution.</returns>
internal static double SampleBeta(Random rnd, double a, double b)
{
double x = Gamma.SampleGamma(rnd, a, 1.0);
double y = Gamma.SampleGamma(rnd, b, 1.0);
return x / (x + y);
}
}
}

504
src/Managed/Distributions/Continuous/Gamma.cs

@ -0,0 +1,504 @@
// <copyright file="Gamma.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://mathnet.opensourcedotnet.info
//
// Copyright (c) 2009 Math.NET
//
// Permission is hereby granted, free of charge, to any person
// obtaining a copy of this software and associated documentation
// files (the "Software"), to deal in the Software without
// restriction, including without limitation the rights to use,
// copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following
// conditions:
//
// The above copyright notice and this permission notice shall be
// included in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
namespace MathNet.Numerics.Distributions
{
using System;
using System.Collections.Generic;
using Properties;
/// <summary>
/// Implements the univariate Gamma distribution. For details about this distribution, see
/// <a href="http://en.wikipedia.org/wiki/Gamma_distribution">Wikipedia - Gamma distribution</a>.
/// </summary>
/// <remarks>
/// <para>The Gamma distribution is parametrized by a shape and inverse scale parameter. When we want
/// to specify a Gamma distribution which is a point distribution we set the shape parameter to be the
/// location of the point distribution and the inverse scale as positive infinity.</para>
/// <para> Random number generation for the Gamma distribution is based on the algorithm in:
/// "A Simple Method for Generating Gamma Variables" - Marsaglia &amp; Tsang
/// ACM Transactions on Mathematical Software, Vol. 26, No. 3, September 2000, Pages 363–372.</para>
/// <para>The distribution will use the <see cref="System.Random"/> by default.
/// Users can get/set the random number generator by using the <see cref="RandomSource"/> property.</para>
/// <para>The statistics classes will check all the incoming parameters whether they are in the allowed
/// range. This might involve heavy computation. Optionally, by setting Control.CheckDistributionParameters
/// to false, all parameter checks can be turned off.</para></remarks>
public class Gamma : IContinuousDistribution
{
/// <summary>
/// Gamma shape parameter.
/// </summary>
private double _shape;
/// <summary>
/// Gamma inverse scale parameter.
/// </summary>
private double _invScale;
/// <summary>
/// Initializes a new instance of the Gamma distribution.
/// </summary>
/// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param>
public Gamma(double shape, double invScale)
{
SetParameters(shape, invScale);
RandomSource = new Random();
}
/// <summary>
/// Constructs a Gamma distribution from a shape and scale parameter. The distribution will
/// be initialized with the default <seealso cref="System.Random"/> random number generator.
/// </summary>
/// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="scale">The scale of the Gamma distribution.</param>
/// <returns>a normal distribution.</returns>
public static Gamma WithShapeScale(double shape, double scale)
{
return new Gamma(shape, 1.0/scale);
}
/// <summary>
/// Constructs a Gamma distribution from a shape and inverse scale parameter. The distribution will
/// be initialized with the default <seealso cref="System.Random"/> random number generator.
/// </summary>
/// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <returns>a normal distribution.</returns>
public static Gamma WithShapeInvScale(double shape, double invScale)
{
return new Gamma(shape, invScale);
}
/// <summary>
/// A string representation of the distribution.
/// </summary>
/// <returns>a string representation of the distribution.</returns>
public override string ToString()
{
return "Gamma(Shape = " + _shape + ", Inverse Scale = " + _invScale + ")";
}
/// <summary>
/// Checks whether the parameters of the distribution are valid.
/// </summary>
/// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <returns>True when the parameters are valid, false otherwise.</returns>
private static bool IsValidParameterSet(double shape, double invScale)
{
if (shape < 0.0 || invScale < 0.0 || Double.IsNaN(shape) || Double.IsNaN(invScale))
{
return false;
}
return true;
}
/// <summary>
/// Sets the parameters of the distribution after checking their validity.
/// </summary>
/// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
_shape = shape;
_invScale = invScale;
}
/// <summary>
/// Gets or sets the shape of the Gamma distribution.
/// </summary>
public double Shape
{
get
{
return _shape;
}
set
{
SetParameters(value, _invScale);
}
}
/// <summary>
/// Gets or sets the scale of the Gamma distribution.
/// </summary>
public double Scale
{
get
{
return 1.0 / _invScale;
}
set
{
SetParameters(_shape, 1.0/value);
}
}
/// <summary>
/// Gets or sets the inverse scale of the Gamma distribution.
/// </summary>
public double InvScale
{
get
{
return _invScale;
}
set
{
SetParameters(_shape, value);
}
}
#region IDistribution implementation
/// <summary>
/// Gets or sets the random number generator which is used to draw random samples.
/// </summary>
public Random RandomSource { get; set; }
/// <summary>
/// Gets the mean of the Gamma distribution.
/// </summary>
public double Mean
{
get
{
if (Double.IsPositiveInfinity(_invScale))
{
return _shape;
}
else
{
return _shape / _invScale;
}
}
}
/// <summary>
/// Gets the variance of the Gamma distribution.
/// </summary>
public double Variance
{
get
{
if (Double.IsPositiveInfinity(_invScale))
{
return 0.0;
}
else
{
return _shape / (_invScale * _invScale);
}
}
}
/// <summary>
/// Gets the standard deviation of the Gamma distribution.
/// </summary>
public double StdDev
{
get
{
if (Double.IsPositiveInfinity(_invScale))
{
return 0.0;
}
else
{
return Math.Sqrt(_shape / (_invScale * _invScale));
}
}
}
/// <summary>
/// Gets the entropy of the Gamma distribution.
/// </summary>
public double Entropy
{
get
{
if (Double.IsPositiveInfinity(_invScale))
{
return 0.0;
}
else
{
return _shape - Math.Log(_invScale) + SpecialFunctions.GammaLn(_shape) + (1.0 - _shape) * SpecialFunctions.DiGamma(_shape);
}
}
}
/// <summary>
/// Gets the skewness of the Gamma distribution.
/// </summary>
public double Skewness
{
get { return 2.0 / Math.Sqrt(_shape); }
}
#endregion
#region IContinuousDistribution implementation
/// <summary>
/// Gets the mode of the Gamma distribution.
/// </summary>
public double Mode
{
get
{
if (Double.IsPositiveInfinity(_invScale))
{
return _shape;
}
else
{
return (_shape - 1.0) / _invScale;
}
}
}
/// <summary>
/// Gets the median of the Gamma distribution.
/// </summary>
public double Median
{
get { throw new NotSupportedException(); }
}
/// <summary>
/// Gets the minimum of the Gamma distribution.
/// </summary>
public double Minimum
{
get { return 0.0; }
}
/// <summary>
/// Gets the maximum of the Gamma distribution.
/// </summary>
public double Maximum
{
get { return Double.PositiveInfinity; }
}
/// <summary>
/// Computes the density of the Gamma distribution.
/// </summary>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
public double Density(double x)
{
if (Double.IsPositiveInfinity(_invScale))
{
if (x == _shape)
{
return Double.PositiveInfinity;
}
else
{
return 0.0;
}
}
else
{
return Math.Pow(_invScale, _shape) * Math.Pow(x, _shape - 1.0) * Math.Exp(-_invScale * x) / SpecialFunctions.Gamma(_shape);
}
}
/// <summary>
/// Computes the log density of the Gamma distribution.
/// </summary>
/// <param name="x">The location at which to compute the log density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
public double DensityLn(double x)
{
if (Double.IsPositiveInfinity(_invScale))
{
if (x == _shape)
{
return Double.PositiveInfinity;
}
else
{
return Double.NegativeInfinity;
}
}
else
{
return _shape * Math.Log(_invScale) + (_shape - 1.0) * Math.Log(x) - _invScale * x - SpecialFunctions.GammaLn(_shape);
}
}
/// <summary>
/// Computes the cumulative distribution function of the Gamma distribution.
/// </summary>
/// <param name="x">The location at which to compute the cumulative density.</param>
/// <returns>the cumulative density at <paramref name="x"/>.</returns>
public double CumulativeDistribution(double x)
{
if (Double.IsPositiveInfinity(_invScale))
{
if (x >= _shape)
{
return 1.0;
}
else
{
return 0.0;
}
}
else
{
return SpecialFunctions.IncompleteGamma(_shape, x * _invScale, true);
}
}
/// <summary>
/// Generates a sample from the Gamma distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return SampleGamma(RandomSource, _shape, _invScale);
}
/// <summary>
/// Generates a sequence of samples from the Gamma distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
while (true)
{
yield return SampleGamma(RandomSource, _shape, _invScale);
}
}
#endregion
/// <summary>
/// Generates a sample from the Gamma distribution.
/// </summary>
/// <param name="rng">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution from which to generate samples.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution from which to generate samples.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rng, double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleGamma(rng, shape, invScale);
}
/// <summary>
/// Generates a sequence of samples from the Gamma distribution.
/// </summary>
/// <param name="rng">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution from which to generate samples.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution from which to generate samples.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rng, double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleGamma(rng, shape, invScale);
}
}
/// <summary>
/// Sampling implementation based on:
/// "A Simple Method for Generating Gamma Variables" - Marsaglia &amp; Tsang
/// ACM Transactions on Mathematical Software, Vol. 26, No. 3, September 2000, Pages 363–372.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <returns>A sample from a Gamma distributed random variable.</returns>
internal static double SampleGamma(System.Random rnd, double shape, double invScale)
{
if (Double.IsPositiveInfinity(invScale))
{
return shape;
}
else
{
double a = shape;
double alphafix = 1.0;
// Fix when alpha is less than one.
if (shape < 1.0)
{
a = shape + 1.0;
alphafix = System.Math.Pow(rnd.NextDouble(), 1.0 / shape);
}
double d = a - 1.0 / 3.0;
double c = 1.0 / System.Math.Sqrt(9.0 * d);
while (true)
{
double x = Normal.Sample(rnd, 0.0, 1.0);
double v = 1.0 + c * x;
while (v <= 0.0)
{
x = Normal.Sample(rnd, 0.0, 1.0);
v = 1.0 + c * x;
}
v = v * v * v;
double u = rnd.NextDouble();
x = x * x;
if (u < 1.0 - 0.0331 * x * x)
{
return alphafix * d * v / invScale;
}
if (System.Math.Log(u) < 0.5 * x + d * (1.0 - v + System.Math.Log(v)))
{
return alphafix * d * v / invScale;
}
}
}
}
}
}

2
src/Managed/Managed.csproj

@ -48,7 +48,9 @@
<Compile Include="Complex.cs" />
<Compile Include="Constants.cs" />
<Compile Include="Control.cs" />
<Compile Include="Distributions\Continuous\Beta.cs" />
<Compile Include="Distributions\Continuous\ContinuousUniform.cs" />
<Compile Include="Distributions\Continuous\Gamma.cs" />
<Compile Include="Distributions\Continuous\Normal.cs" />
<Compile Include="Distributions\Discrete\Bernoulli.cs" />
<Compile Include="Distributions\IContinuousDistribution.cs" />

32
src/Managed/SpecialFunctions.cs

@ -58,5 +58,37 @@ namespace MathNet.Numerics
return 0d;
}
public static double BetaLn(double a, double b)
{
return Double.NaN;
}
public static double BetaRegularized(double a, double b, double x)
{
return Double.NaN;
}
public static double DiGamma(double x)
{
return Double.NaN;
}
public static double Gamma(double x)
{
return Double.NaN;
}
public static double GammaLn(double x)
{
return Double.NaN;
}
public static double IncompleteGamma(double x, double z, bool reg)
{
return Double.NaN;
}
}
}

6
src/Native.UnitTests/Native.UnitTests.csproj

@ -68,9 +68,15 @@
<Compile Include="..\Managed.UnitTests\ComplexTests\ComplexTest.cs">
<Link>ComplexTests\ComplexTest.cs</Link>
</Compile>
<Compile Include="..\Managed.UnitTests\DistributionTests\Continuous\BetaTests.cs">
<Link>DistributionTests\Continuous\BetaTests.cs</Link>
</Compile>
<Compile Include="..\Managed.UnitTests\DistributionTests\Continuous\ContinuousUniformTests.cs">
<Link>DistributionTests\Continuous\ContinuousUniformTests.cs</Link>
</Compile>
<Compile Include="..\Managed.UnitTests\DistributionTests\Continuous\GammaTests.cs">
<Link>DistributionTests\Continuous\GammaTests.cs</Link>
</Compile>
<Compile Include="..\Managed.UnitTests\DistributionTests\Continuous\NormalTests.cs">
<Link>DistributionTests\Continuous\NormalTests.cs</Link>
</Compile>

6
src/Native/Native.csproj

@ -53,9 +53,15 @@
<Compile Include="..\Managed\Control.cs">
<Link>Control.cs</Link>
</Compile>
<Compile Include="..\Managed\Distributions\Continuous\Beta.cs">
<Link>Distributions\Continuous\Beta.cs</Link>
</Compile>
<Compile Include="..\Managed\Distributions\Continuous\ContinuousUniform.cs">
<Link>Distributions\Continuous\ContinuousUniform.cs</Link>
</Compile>
<Compile Include="..\Managed\Distributions\Continuous\Gamma.cs">
<Link>Distributions\Continuous\Gamma.cs</Link>
</Compile>
<Compile Include="..\Managed\Distributions\Continuous\Normal.cs">
<Link>Distributions\Continuous\Normal.cs</Link>
</Compile>

Loading…
Cancel
Save