Migrated inverse gamma distribution to shape fitter

This commit is contained in:
Thorsten Sommer 2020-10-31 15:43:08 +01:00
parent 8b752c17ee
commit 66d58e74d3
3 changed files with 74 additions and 70 deletions

View File

@ -1,36 +0,0 @@
using System;
using System.Threading;
using System.Threading.Tasks;
namespace FastRng.Double.Distributions
{
public sealed class InverseGamma : IDistribution
{
private double shape = 1.0;
public IRandom Random { get; set; }
public double Shape
{
get => this.shape;
set
{
if(value <= 0.0)
throw new ArgumentOutOfRangeException(message: "Shape must be greater than 0", null);
this.shape = value;
}
}
public double Scale { get; set; } = 1.0;
public async ValueTask<double> GetDistributedValue(CancellationToken token = default)
{
if (this.Random == null)
return double.NaN;
var gammaDist = new Gamma{ Shape = this.Shape, Scale = 1.0 / this.Scale };
return 1.0 / await this.Random.NextNumber(gammaDist, token);
}
}
}

View File

@ -0,0 +1,46 @@
using System;
using System.Threading;
using System.Threading.Tasks;
namespace FastRng.Double.Distributions
{
public sealed class InverseGammaA3B05 : IDistribution
{
private const double ALPHA = 3.0;
private const double BETA = 0.5;
private const double CONSTANT = 0.213922656884911;
private static readonly double FACTOR_LEFT;
private ShapeFitter fitter;
private IRandom random;
static InverseGammaA3B05()
{
var gammaAlpha = MathTools.Gamma(ALPHA);
var betaToTheAlpha = Math.Pow(BETA, ALPHA);
FACTOR_LEFT = CONSTANT * (betaToTheAlpha / gammaAlpha);
}
public IRandom Random
{
get => this.random;
set
{
this.random = value;
this.fitter = new ShapeFitter(InverseGammaA3B05.ShapeFunction, this.random, 100);
}
}
private static double ShapeFunction(double x) => FACTOR_LEFT * Math.Pow(x, -ALPHA - 1.0d) * Math.Exp(-BETA / x);
public async ValueTask<double> GetDistributedValue(CancellationToken token = default)
{
if (this.Random == null)
return double.NaN;
return await this.fitter.NextNumber(token);
}
}
}

View File

@ -8,31 +8,40 @@ using NUnit.Framework;
namespace FastRngTests.Double.Distributions namespace FastRngTests.Double.Distributions
{ {
[ExcludeFromCodeCoverage] [ExcludeFromCodeCoverage]
public class InverseGamma public class InverseGammaA3B05
{ {
[Test] [Test]
[Category(TestCategories.COVER)] [Category(TestCategories.COVER)]
[Category(TestCategories.NORMAL)] [Category(TestCategories.NORMAL)]
public async Task TestInverseGammaDistribution01() public async Task TestInverseGammaDistribution01()
{ {
const double SHAPE = 3.0; var dist = new FastRng.Double.Distributions.InverseGammaA3B05();
const double SCALE = 1.2; var fra = new FrequencyAnalysis();
const double MEAN = SCALE / (SHAPE - 1);
const double VARIANCE = SCALE * SCALE / ((SHAPE - 1) * (SHAPE - 1) * (SHAPE - 2));
var dist = new FastRng.Double.Distributions.InverseGamma{ Shape = SHAPE, Scale = SCALE };
var stats = new RunningStatistics();
var rng = new MultiThreadedRng(); var rng = new MultiThreadedRng();
for (var n = 0; n < 100_000; n++) for (var n = 0; n < 100_000; n++)
stats.Push(await rng.NextNumber(dist)); fra.CountThis(await rng.NextNumber(dist));
rng.StopProducer(); rng.StopProducer();
TestContext.WriteLine($"mean={MEAN} vs. {stats.Mean}"); var result = fra.NormalizeAndPlotEvents(TestContext.WriteLine);
TestContext.WriteLine($"variance={VARIANCE} vs {stats.Variance}");
Assert.That(result[0], Is.EqualTo(0.0000000000000003).Within(0.0000001));
Assert.That(result[1], Is.EqualTo(0.0000011605257228).Within(0.00001));
Assert.That(result[2], Is.EqualTo(0.0009536970016103).Within(0.0005));
Assert.That(stats.Mean, Is.EqualTo(MEAN).Within(0.1), "Mean is out of range"); Assert.That(result[21], Is.EqualTo(0.5880485243048120).Within(0.05));
Assert.That(stats.Variance, Is.EqualTo(VARIANCE).Within(0.1), "Variance is out of range"); Assert.That(result[22], Is.EqualTo(0.5433842148912880).Within(0.05));
Assert.That(result[23], Is.EqualTo(0.5017780549216030).Within(0.05));
Assert.That(result[50], Is.EqualTo(0.0741442015957425).Within(0.006));
Assert.That(result[75], Is.EqualTo(0.0207568945092484).Within(0.006));
Assert.That(result[85], Is.EqualTo(0.0136661506653688).Within(0.004));
Assert.That(result[90], Is.EqualTo(0.0112550619601327).Within(0.004));
Assert.That(result[97], Is.EqualTo(0.0087026933539773).Within(0.002));
Assert.That(result[98], Is.EqualTo(0.0083995375385004).Within(0.002));
Assert.That(result[99], Is.EqualTo(0.0081094156379928).Within(0.002));
} }
[Test] [Test]
@ -41,9 +50,10 @@ namespace FastRngTests.Double.Distributions
public async Task TestInverseGammaGeneratorWithRange01() public async Task TestInverseGammaGeneratorWithRange01()
{ {
var rng = new MultiThreadedRng(); var rng = new MultiThreadedRng();
var dist = new FastRng.Double.Distributions.InverseGammaA3B05();
var samples = new double[1_000]; var samples = new double[1_000];
for (var n = 0; n < samples.Length; n++) for (var n = 0; n < samples.Length; n++)
samples[n] = await rng.NextNumber(-1.0, 1.0, new FastRng.Double.Distributions.InverseGamma()); samples[n] = await rng.NextNumber(-1.0, 1.0, dist);
rng.StopProducer(); rng.StopProducer();
Assert.That(samples.Min(), Is.GreaterThanOrEqualTo(-1.0), "Min out of range"); Assert.That(samples.Min(), Is.GreaterThanOrEqualTo(-1.0), "Min out of range");
@ -56,9 +66,10 @@ namespace FastRngTests.Double.Distributions
public async Task TestInverseGammaGeneratorWithRange02() public async Task TestInverseGammaGeneratorWithRange02()
{ {
var rng = new MultiThreadedRng(); var rng = new MultiThreadedRng();
var dist = new FastRng.Double.Distributions.InverseGammaA3B05();
var samples = new double[1_000]; var samples = new double[1_000];
for (var n = 0; n < samples.Length; n++) for (var n = 0; n < samples.Length; n++)
samples[n] = await rng.NextNumber(0.0, 1.0, new FastRng.Double.Distributions.InverseGamma()); samples[n] = await rng.NextNumber(0.0, 1.0, dist);
rng.StopProducer(); rng.StopProducer();
Assert.That(samples.Min(), Is.GreaterThanOrEqualTo(0.0), "Min is out of range"); Assert.That(samples.Min(), Is.GreaterThanOrEqualTo(0.0), "Min is out of range");
@ -71,7 +82,7 @@ namespace FastRngTests.Double.Distributions
public async Task TestInverseGammaGeneratorWithRange03() public async Task TestInverseGammaGeneratorWithRange03()
{ {
var rng = new MultiThreadedRng(); var rng = new MultiThreadedRng();
var dist = new FastRng.Double.Distributions.InverseGamma { Random = rng }; // Test default parameters var dist = new FastRng.Double.Distributions.InverseGammaA3B05 { Random = rng }; // Test default parameters
var samples = new double[1_000]; var samples = new double[1_000];
for (var n = 0; n < samples.Length; n++) for (var n = 0; n < samples.Length; n++)
@ -82,29 +93,12 @@ namespace FastRngTests.Double.Distributions
Assert.That(samples.Max(), Is.LessThanOrEqualTo(1.0), "Max is out of range"); Assert.That(samples.Max(), Is.LessThanOrEqualTo(1.0), "Max is out of range");
} }
[Test]
[Category(TestCategories.COVER)]
[Category(TestCategories.NORMAL)]
public void ParameterTest01()
{
var dist = new FastRng.Double.Distributions.InverseGamma();
Assert.Throws<ArgumentOutOfRangeException>(() => dist.Shape = 0);
Assert.Throws<ArgumentOutOfRangeException>(() => dist.Shape = -78);
Assert.DoesNotThrow(() => dist.Shape = 0.0001);
Assert.DoesNotThrow(() => dist.Shape = 4);
Assert.DoesNotThrow(() => dist.Scale = -45);
Assert.DoesNotThrow(() => dist.Scale = 15);
Assert.DoesNotThrow(() => dist.Scale = 0);
}
[Test] [Test]
[Category(TestCategories.COVER)] [Category(TestCategories.COVER)]
[Category(TestCategories.NORMAL)] [Category(TestCategories.NORMAL)]
public async Task NoRandomNumberGenerator01() public async Task NoRandomNumberGenerator01()
{ {
var dist = new FastRng.Double.Distributions.InverseGamma(); var dist = new FastRng.Double.Distributions.InverseGammaA3B05();
Assert.DoesNotThrowAsync(async () => await dist.GetDistributedValue()); Assert.DoesNotThrowAsync(async () => await dist.GetDistributedValue());
Assert.That(await dist.GetDistributedValue(), Is.NaN); Assert.That(await dist.GetDistributedValue(), Is.NaN);
} }