Browse Source

Added new distributions and corresponding tests

- New distributions are
   - Three parameter Burr (Type XII)
   - Inverse Gaussian
   - Truncated Pareto
 - Unit tests added for all public methods
 - New distributions added to CommonDistributionTests
uap-experimental
Philip Gruber 7 years ago
parent
commit
9028337af9
  1. 3
      src/Numerics.Tests/DistributionTests/CommonDistributionTests.cs
  2. 426
      src/Numerics.Tests/DistributionTests/Continuous/BurrTests.cs
  3. 419
      src/Numerics.Tests/DistributionTests/Continuous/InverseGaussianTests.cs
  4. 446
      src/Numerics.Tests/DistributionTests/Continuous/TruncatedParetoTests.cs
  5. 429
      src/Numerics/Distributions/Burr.cs
  6. 449
      src/Numerics/Distributions/InverseGaussian.cs
  7. 466
      src/Numerics/Distributions/TruncatedPareto.cs
  8. 9
      src/Numerics/Precision.Comparison.cs
  9. 11
      src/Numerics/Properties/Resources.Designer.cs
  10. 3
      src/Numerics/Properties/Resources.resx

3
src/Numerics.Tests/DistributionTests/CommonDistributionTests.cs

@ -67,6 +67,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
{ {
new Beta(1.0, 1.0), new Beta(1.0, 1.0),
new BetaScaled(1.0, 1.5, 0.5, 2.0), new BetaScaled(1.0, 1.5, 0.5, 2.0),
new Burr(2.0, 3.0, 1.0),
new Cauchy(1.0, 1.0), new Cauchy(1.0, 1.0),
new Chi(3.0), new Chi(3.0),
new ChiSquared(3.0), new ChiSquared(3.0),
@ -76,6 +77,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
new FisherSnedecor(0.3, 0.4), new FisherSnedecor(0.3, 0.4),
new Gamma(1.0, 1.0), new Gamma(1.0, 1.0),
new InverseGamma(1.0, 1.0), new InverseGamma(1.0, 1.0),
new InverseGaussian(1.0, 3.0),
new Laplace(1.0, 0.5), new Laplace(1.0, 0.5),
new LogNormal(1.0, 1.0), new LogNormal(1.0, 1.0),
new Normal(0.0, 1.0), new Normal(0.0, 1.0),
@ -84,6 +86,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
new Stable(0.5, 1.0, 0.5, 1.0), new Stable(0.5, 1.0, 0.5, 1.0),
new StudentT(0.0, 1.0, 5.0), new StudentT(0.0, 1.0, 5.0),
new Triangular(0, 1, 0.7), new Triangular(0, 1, 0.7),
new TruncatedPareto(1.0, 3.0, 1000),
new Weibull(1.0, 1.0), new Weibull(1.0, 1.0),
}; };

426
src/Numerics.Tests/DistributionTests/Continuous/BurrTests.cs

@ -0,0 +1,426 @@
// <copyright file="BurrTests.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
//
// Copyright (c) 2009-2019 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>
using MathNet.Numerics.Distributions;
using NUnit.Framework;
using System;
using System.Linq;
namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
{
/// <summary>
/// <c>Burr</c> distribution tests.
/// </summary>
[TestFixture, Category("Distributions")]
public class BurrTests
{
private const int highPrecision = 12;
private const int lowPrecision = 8;
/// <summary>
/// Can create <c>Burr</c>.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
[TestCase(1.0, 1.0, 0.5)]
[TestCase(10.0, 0.1, 0.2)]
[TestCase(5.0, 1.0, 3)]
[TestCase(2.0, 10.0, 2.0)]
[TestCase(10.0, 100.0, 1.0)]
public void CanCreateBurr(double a, double c, double k)
{
var n = new Burr(a, c, k);
Assert.AreEqual(a, n.a);
Assert.AreEqual(c, n.c);
Assert.AreEqual(k, n.k);
}
/// <summary>
/// Can create <c>Burr</c> with random source.
/// </summary>
[Test]
public void CanCreateBurrWithRandomSource()
{
var randomSource = new Numerics.Random.MersenneTwister(100);
var n = new Burr(1.0, 1.0, 0.5, randomSource);
Assert.AreEqual(randomSource, n.RandomSource);
}
/// <summary>
/// <c>Burr</c> create fails with bad parameters.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
[TestCase(5.0, 1.0, -1.0)]
[TestCase(-1.0, 1.0, 2.05)]
[TestCase(1.0, -1.0, 2.05)]
[TestCase(Double.NaN, Double.NaN, Double.NaN)]
[TestCase(Double.NaN, 1.0, 2.4)]
[TestCase(1.0, 1.0, double.PositiveInfinity)]
[TestCase(1.0, -1.0, Double.NegativeInfinity)]
public void BurrCreateFailsWithBadParameters(double a, double c, double k)
{
Assert.That(() => new Burr(a, c, k), Throws.ArgumentException);
}
/// <summary>
/// Validate IsValidParameterSet.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="validity">exoected validity of paremeter set.</param>
[TestCase(1.0, 1.0, 0.5, true)]
[TestCase(10.0, 0.1, 0.2, true)]
[TestCase(5.0, 1.0, -1.0, false)]
[TestCase(1.0, 1.0, double.PositiveInfinity, false)]
public void ValidateIsValidParameterSet(double a, double c, double k, bool validity)
{
Assert.AreEqual(Burr.IsValidParameterSet(a, c, k), validity);
}
/// <summary>
/// Validate ToString.
/// </summary>
[Test]
public void ValidateToString()
{
var n = new Burr(1d, 2d, 3d);
Assert.AreEqual("Burr(a = 1, c = 2, k = 3)", n.ToString());
}
/// <summary>
/// Validate mode.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="mode">Expected value.</param>
[TestCase(2.000000, 2.00000, 1.0, 1.154700538379251)]
[TestCase(1.000000, 1.00000, 0.5, 0.0)]
[TestCase(3.000000, 4.00000, 1.0, 2.640335210380180)]
public void ValidateMode(double a, double c, double k, double mode)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(mode, n.Mode, highPrecision);
}
/// <summary>
/// Validate median.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="median">Expected value.</param>
[TestCase(1.000000, 1.00000, 0.5, 3.0)]
[TestCase(1.000000, 0.50000, 0.5, 9.0)]
[TestCase(1.000000, 1.00000, 5.0, 0.148698354997035)]
public void ValidateMedian(double a, double c, double k, double median)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(median, n.Median, highPrecision);
}
/// <summary>
/// Validate mean.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="mean">Expected value.</param>
[TestCase(1.000000, 1.00000, 1.5, 2)]
[TestCase(4.000000, 5.00000, 0.5, 6.198785110989412)]
[TestCase(4.000000, 5.00000, 5.0, 2.729694550490384)]
public void ValidateMean(double a, double c, double k, double mean)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(mean, n.Mean, highPrecision);
}
/// <summary>
/// Validate variance.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="variance">Expected value.</param>
[TestCase(4.0, 5.0, 2.0, 0.983559126843161)]
[TestCase(2.0, 3.5, 4.0, 0.207489170174404)]
public void ValidateVariance(double a, double c, double k, double variance)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(variance, n.Variance, highPrecision);
}
/// <summary>
/// Validate standard deviation.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="std">Expected value.</param>
[TestCase(4.0, 5.0, 2.0, 0.991745494995143)]
[TestCase(2.0, 3.5, 4.0, 0.455509791524182)]
public void ValidateStandardDeviation(double a, double c, double k, double std)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(std, n.StdDev, highPrecision);
}
/// <summary>
/// Validate minimum.
/// </summary>
[Test]
public void ValidateMinimum()
{
var n = new Burr(1.0, 2.0, 2.0);
Assert.AreEqual(0.0, n.Minimum);
}
/// <summary>
/// Validate maximum.
/// </summary>
[Test]
public void ValidateMaximum()
{
var n = new Burr(1.0, 2.0, 1.0);
Assert.AreEqual(Double.PositiveInfinity, n.Maximum);
}
/// <summary>
/// Validate entropy throws not supported exception
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
[TestCase(1.0, 1.0, 0.5)]
public void ValidateEntropyFailsWithNotSupported(double a, double c, double k)
{
var n = new Burr(a, c, k);
Assert.Throws<NotSupportedException>(() => { var x = n.Entropy; });
}
/// <summary>
/// Validate GetMoments.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="n"> order of the moment.</param>
/// <param name="value">Expected value.</param>
[TestCase(1.000000, 1.00000, 1.5, 1, 2)]
[TestCase(4.000000, 5.00000, 0.5, 1, 6.198785110989412)]
[TestCase(4.000000, 5.00000, 5.0, 1, 2.729694550490384)]
[TestCase(4.0, 5.0, 2.0, 3, 50.738165747621750)]
[TestCase(2.0, 3.5, 4.0, 3, 2.895046685294824)]
public void ValidateMoments(double a, double c, double k, int order, double value)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(value, n.GetMoment(order), lowPrecision);
}
/// <summary>
/// Validate GetMoments throws when given bad parameters
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="order">the order of the moment.</param>
[TestCase(1.0, 1.0, 0.5, 1)]
[TestCase(10.0, 1.0, 1.5, 2)]
[TestCase(5.0, 2, 1.4, 3)]
public void ValidateGetMomentsFailsWithBadParameters(double a, double c, double k, double order)
{
var n = new Burr(a, c, k);
Assert.That(() => n.GetMoment(order), Throws.ArgumentException);
}
/// <summary>
/// Validate skewness.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="skewness">Expected value.</param>
[TestCase(4.0, 5.0, 2.0, 0.635277200891842)]
[TestCase(2.0, 3.5, 4.0, 0.483073007212360)]
public void ValidateSkewness(double a, double c, double k, double skewness)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(skewness, n.Skewness, lowPrecision);
}
/// <summary>
/// Validate density.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="x">Input X value.</param>
/// <param name="p">Expected value.</param>
[TestCase(1.00000, 1.000000, 0.5, 0.100000, 0.433392086020724)]
[TestCase(1.000000, 1.000000, 5.0, 0.500000, 0.438957475994512)]
public void ValidateDensity(double a, double c, double k, double x, double p)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(p, n.Density(x), highPrecision);
AssertHelpers.AlmostEqualRelative(p, Burr.PDF(a, c, k, x), highPrecision);
}
/// <summary>
/// Validate density log.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="x">Input X value.</param>
/// <param name="p">Expected value.</param>
[TestCase(1.000000, 1.000000, 0.5, 0.500000, -1.301344842722192)]
[TestCase(1.000000, 1.000000, 5.0, 0.500000, -0.823352736214886)]
[TestCase(4.000000, 5.000000, 1.0, 5.000000, -1.682583872895781)]
public void ValidateDensityLn(double a, double c, double k, double x, double p)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(p, n.DensityLn(x), highPrecision);
AssertHelpers.AlmostEqualRelative(p, Burr.PDFLn(a, c, k, x), highPrecision);
}
/// <summary>
/// Validate cumulative distribution.
/// </summary>
/// <param name="a">a parameter.</param>
/// <param name="c">c parameter.</param>
/// <param name="k">k parameter.</param>
/// <param name="x">The location at which to compute the density.</param>
/// <param name="f">Expected value.</param>
[TestCase(1.000000, 1.000000, 0.5, 0.500000, 0.183503419072274)]
[TestCase(4.000000, 5.000000, 0.5, 0.500000, 0.00001525843981653452)]
[TestCase(4.000000, 5.000000, 0.5, 5.00000, 0.503203804978536)]
[TestCase(4.000000, 5.000000, 5.0, 5.00000, 0.999084238019496)]
public void ValidateCumulativeDistribution(double a, double c, double k, double x, double f)
{
var n = new Burr(a, c, k);
AssertHelpers.AlmostEqualRelative(f, n.CumulativeDistribution(x), lowPrecision);
AssertHelpers.AlmostEqualRelative(f, Burr.CDF(a, c, k, x), lowPrecision);
}
/// <summary>
/// Can sample static.
/// </summary>
[Test]
public void CanSampleStatic()
{
Burr.Sample(new Numerics.Random.MersenneTwister(100), 1.0, 2.0, 2.0);
}
/// <summary>
/// Can fill array with samples static.
/// </summary>
[Test]
public void CanFillSampleArrayStatic()
{
double[] samples = new double[100];
Burr.Samples(new Numerics.Random.MersenneTwister(100), samples, 1.0, 2.0, 2.0);
Assert.IsTrue(!samples.Any(x => x == 0));
}
/// <summary>
/// Can sample sequence static.
/// </summary>
[Test]
public void CanSampleSequenceStatic()
{
var ied = Burr.Samples(new Numerics.Random.MersenneTwister(100), 1.0, 2.0, 2.0);
GC.KeepAlive(ied.Take(5).ToArray());
}
/// <summary>
/// Fail sample static with bad parameters.
/// </summary>
[Test]
public void FailSampleStatic()
{
Assert.That(() => { var d = Burr.Sample(new Numerics.Random.MersenneTwister(100), 1.0, -1.0, 2.0); }, Throws.ArgumentException);
}
/// <summary>
/// Fail filling array with samples with bad parameters static.
/// </summary>
[Test]
public void FailFillingSampleArrayStatic()
{
double[] samples = new double[100];
Assert.That(() => { Burr.Samples(new Numerics.Random.MersenneTwister(100), samples, -1.0, 2.0, 2.0); }, Throws.ArgumentException);
}
/// <summary>
/// Fail sample sequence static with bad parameters.
/// </summary>
[Test]
public void FailSampleSequenceStatic()
{
Assert.That(() => { var ied = Burr.Samples(new Numerics.Random.MersenneTwister(100), 1.0, -1.0, 2.0).First(); }, Throws.ArgumentException);
}
/// <summary>
/// Can sample.
/// </summary>
[Test]
public void CanSample()
{
var n = new Burr(1.0, 2.0, 2.0);
n.Sample();
}
/// <summary>
/// Can fill array with samples.
/// </summary>
[Test]
public void CanFillSampleArray()
{
double[] samples = new double[100];
var n = new Burr(1.0, 2.0, 2.0, new Numerics.Random.MersenneTwister(100));
n.Samples(samples);
Assert.IsTrue(!samples.Any(x => x == 0));
}
/// <summary>
/// Can sample sequence.
/// </summary>
[Test]
public void CanSampleSequence()
{
var n = new Burr(1.0, 2.0, 2.0);
var ied = n.Samples();
GC.KeepAlive(ied.Take(5).ToArray());
}
}
}

419
src/Numerics.Tests/DistributionTests/Continuous/InverseGaussianTests.cs

@ -0,0 +1,419 @@
// <copyright file="InverseGaussianTests.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
//
// Copyright (c) 2009-2019 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>
using MathNet.Numerics.Distributions;
using NUnit.Framework;
using System;
using System.Linq;
namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
{
/// <summary>
/// <c>InverseGaussian</c> distribution tests.
/// </summary>
[TestFixture, Category("Distributions")]
public class InverseGaussianTests
{
private const int precision = 12;
/// <summary>
/// Can create <c>Inverse Gaussian</c>.
/// </summary>
/// <param name="mu">Mu value.</param>
/// <param name="lambda">Lambda parameter.</param>
[TestCase(0.1, 0.1)]
[TestCase(1.0, 1.0)]
[TestCase(10.0, 10.0)]
public void CanCreateInverseGaussian(double mu, double lambda)
{
var n = new InverseGaussian(mu, lambda);
Assert.AreEqual(mu, n.Mu);
Assert.AreEqual(lambda, n.Lambda);
}
/// <summary>
/// Can create <c>Inverse Gaussian</c> with random source.
/// </summary>
[Test]
public void CanCreateInverseGaussianWithRandomSource()
{
var randomSource = new Numerics.Random.MersenneTwister(100);
var n = new InverseGaussian(1.0, 1.0, randomSource);
Assert.AreEqual(randomSource, n.RandomSource);
}
/// <summary>
/// <c>InverseGaussian</c> create fails with bad parameters.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
[TestCase(Double.NaN, 1.0)]
[TestCase(1.0, Double.NaN)]
[TestCase(Double.NaN, Double.NaN)]
[TestCase(-1.0, -1.0)]
public void InverseGaussianCreateFailsWithBadParameters(double mu, double lambda)
{
Assert.That(() => new LogNormal(mu, lambda), Throws.ArgumentException);
}
/// <summary>
/// Validate IsValidParameterSet.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
/// <param name="validity">exoected validity of paremeter set.</param>
[TestCase(1.0, 1.0, true)]
[TestCase(10.0, 10.0, true)]
[TestCase(-1.0, -1.0, false)]
public void ValidateIsValidParameterSet(double mu, double lambda, bool validity)
{
Assert.AreEqual(InverseGaussian.IsValidParameterSet(mu, lambda), validity);
}
/// <summary>
/// Validate ToString.
/// </summary>
[Test]
public void ValidateToString()
{
var n = new InverseGaussian(1d, 2d);
Assert.AreEqual("InverseGaussian(μ = 1, λ = 2)", n.ToString());
}
/// <summary>
/// Validate mode.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
/// <param name="mode">Expected value.</param>
[TestCase(0.100000, 0.100000, 0.030277563773199)]
[TestCase(1.500000, 5.500000, 1.007026953273370)]
[TestCase(2.500000, 1.500000, 0.481456008918130)]
[TestCase(5.500000, 2.500000, 0.815033614523337)]
[TestCase(5.500000, 5.500000, 1.665266007525970)]
public void ValidateMode(double mu, double lambda, double mode)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(mode, n.Mode, precision);
}
/// <summary>
/// Validate median.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
/// <param name="median">Expected value.</param>
[TestCase(0.100000, 0.100000, 0.067584130569524)]
[TestCase(1.500000, 0.100000, 0.190369896722453)]
[TestCase(2.500000, 0.100000, 0.201110669500521)]
[TestCase(5.500000, 5.500000, 3.717127181323815)]
public void ValidateMedian(double mu, double lambda, double median)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(median, n.Median, precision);
}
/// <summary>
/// Validate mean.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
/// <param name="mean">Expected value.</param>
[TestCase(0.100000, 5.500000, 0.1)]
[TestCase(1.500000, 0.100000, 1.5)]
[TestCase(2.500000, 5.500000, 2.5)]
[TestCase(5.500000, 0.100000, 5.5)]
public void ValidateMean(double mu, double lambda, double mean)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(mean, n.Mean, precision);
}
/// <summary>
/// Validate variance.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
/// <param name="variance">Expected value.</param>
[TestCase(0.100000, 5.500000, 1.818181818181819e-04)]
[TestCase(1.500000, 0.100000, 33.75)]
public void ValidateVariance(double mu, double lambda, double variance)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(variance, n.Variance, precision);
}
/// <summary>
/// Validate standard deviation.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
/// <param name="std">Expected value.</param>
[TestCase(2.500000, 5.500000, 1.685499656158105)]
[TestCase(5.500000, 0.100000, 40.789091679026143)]
public void ValidateStandardDeviation(double mu, double lambda, double std)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(std, n.StdDev, precision);
}
/// <summary>
/// Validate minimum.
/// </summary>
[Test]
public void ValidateMinimum()
{
var n = new InverseGaussian(1.0, 2.0);
Assert.AreEqual(0.0, n.Minimum);
}
/// <summary>
/// Validate maximum.
/// </summary>
[Test]
public void ValidateMaximum()
{
var n = new InverseGaussian(1.0, 2.0);
Assert.AreEqual(Double.PositiveInfinity, n.Maximum);
}
/// <summary>
/// Validate entropy throws not supported exception
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
[TestCase(1.0, 1.0)]
public void ValidateEntropyFailsWithNotSupported(double mu, double lambda)
{
var n = new InverseGaussian(mu, lambda);
Assert.Throws<NotSupportedException>(() => { var x = n.Entropy; });
}
/// <summary>
/// Validate skewness.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="´lambda">Lambda parameter.</param>
/// <param name="skewness">Expected value.</param>
[TestCase(0.100000, 1.500000, 0.774596669241483)]
[TestCase(1.500000, 2.500000, 2.323790007724450)]
[TestCase(1.500000, 5.500000, 1.566698903601281)]
[TestCase(2.500000, 0.100000, 15.0)]
[TestCase(2.500000, 1.500000, 3.872983346207417)]
[TestCase(5.500000, 2.500000, 4.449719092257398)]
[TestCase(5.500000, 5.500000, 3.0)]
public void ValidateSkewness(double mu, double lambda, double skewness)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(skewness, n.Skewness, precision);
}
/// <summary>
/// Validate density.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
/// <param name="x">The location at which to compute the density.</param>
/// <param name="p">Expected value.</param>
[TestCase(1.500000, 2.500000, 0.500000, 0.587321148416470)]
[TestCase(1.500000, 2.500000, 0.800000, 0.627284170789435)]
[TestCase(2.500000, 2.500000, 0.500000, 0.360208446721537)]
[TestCase(2.500000, 2.500000, 0.800000, 0.428023216678204)]
public void ValidateDensity(double mu, double lambda, double x, double p)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(p, n.Density(x), precision);
AssertHelpers.AlmostEqualRelative(p, InverseGaussian.PDF(mu, lambda, x), precision);
}
/// <summary>
/// Validate density log.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Sigma value.</param>
/// <param name="x">The location at which to compute the log density.</param>
/// <param name="p">Expected value.</param>
[TestCase(1.500000, 0.100000, 0.100000, 0.948091004233817)]
[TestCase(1.500000, 1.500000, 0.800000, -0.585657318845943)]
[TestCase(2.500000, 0.100000, 0.800000, -1.764415752730381)]
[TestCase(2.500000, 1.500000, 0.100000, -4.174328339659523)]
[TestCase(2.500000, 1.500000, 0.500000, -0.636485208310673)]
public void ValidateDensityLn(double mu, double lambda, double x, double p)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(p, n.DensityLn(x), precision);
AssertHelpers.AlmostEqualRelative(p, InverseGaussian.PDFLn(mu, lambda, x), precision);
}
/// <summary>
/// Validate cumulative distribution.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Sigma value.</param>
/// <param name="x">The location at which to compute the density.</param>
/// <param name="f">Expected value.</param>
[TestCase(1.500000, 2.500000, 0.100000, 0.000002882116410867029)]
[TestCase(2.500000, 1.500000, 0.100000, 0.0001938001952257318)]
[TestCase(2.500000, 1.500000, 0.500000, 0.145457623130791)]
[TestCase(2.500000, 2.500000, 0.100000, 0.000001529605202470422)]
[TestCase(2.500000, 2.500000, 0.800000, 0.187168922781367)]
public void ValidateCumulativeDistribution(double mu, double lambda, double x, double f)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(f, n.CumulativeDistribution(x), precision);
AssertHelpers.AlmostEqualRelative(f, InverseGaussian.CDF(mu, lambda, x), precision);
}
/// <summary>
/// Validate inverse cumulative distribution.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Sigma value.</param>
/// <param name="probability">The location at which to compute the inverse cumulative distribution function.</param>
/// <param name="f">Expected value.</param>
[TestCase(1.500000, 1.500000, 0.100000, 0.356437063090717)]
[TestCase(1.500000, 1.500000, 0.500000, 1.013761958542859)]
[TestCase(2.500000, 2.500000, 0.500000, 1.689603264238098)]
[TestCase(2.500000, 2.500000, 0.800000, 3.619719792074686)]
public void ValidateInverseCumulativeDistribution(double mu, double lambda, double probability, double f)
{
var n = new InverseGaussian(mu, lambda);
AssertHelpers.AlmostEqualRelative(f, InverseGaussian.ICDF(mu, lambda, probability), precision);
AssertHelpers.AlmostEqualRelative(f, n.InvCDF(probability), precision);
}
/// <summary>
/// Can sample static.
/// </summary>
[Test]
public void CanSampleStatic()
{
InverseGaussian.Sample(new Numerics.Random.MersenneTwister(100), 1.0, 1.0);
}
/// <summary>
/// Can fill array with samples static.
/// </summary>
[Test]
public void CanFillSampleArrayStatic()
{
double[] samples = new double[100];
InverseGaussian.Samples(new Numerics.Random.MersenneTwister(100), samples, 1.0, 1.0);
Assert.IsTrue(!samples.Any(x => x == 0));
}
/// <summary>
/// Can sample sequence static.
/// </summary>
[Test]
public void CanSampleSequenceStatic()
{
var ied = InverseGaussian.Samples(new Numerics.Random.MersenneTwister(100), 1.0, 1.0);
GC.KeepAlive(ied.Take(5).ToArray());
}
/// <summary>
/// Fail sample static with bad parameters.
/// </summary>
[Test]
public void FailSampleStatic()
{
Assert.That(() => { var d = InverseGaussian.Sample(new Numerics.Random.MersenneTwister(100), 0.0, -1.0); }, Throws.ArgumentException);
}
/// <summary>
/// Fail filling array with samples with bad parameters static.
/// </summary>
[Test]
public void FailFillingSampleArrayStatic()
{
double[] samples = new double[100];
Assert.That(() => { InverseGaussian.Samples(new Numerics.Random.MersenneTwister(100), samples, -1.0, 1.0); }, Throws.ArgumentException);
}
/// <summary>
/// Fail sample sequence static with bad parameters.
/// </summary>
[Test]
public void FailSampleSequenceStatic()
{
Assert.That(() => { var ied = InverseGaussian.Samples(new Numerics.Random.MersenneTwister(100), 0.0, -1.0).First(); }, Throws.ArgumentException);
}
/// <summary>
/// Can sample.
/// </summary>
[Test]
public void CanSample()
{
var n = new InverseGaussian(1.0, 2.0);
n.Sample();
}
/// <summary>
/// Can fill array with samples.
/// </summary>
[Test]
public void CanFillSampleArray()
{
double[] samples = new double[100];
var n = new InverseGaussian(1.0, 2.0, new Numerics.Random.MersenneTwister(100));
n.Samples(samples);
Assert.IsTrue(!samples.Any(x => x == 0));
}
/// <summary>
/// Can sample sequence.
/// </summary>
[Test]
public void CanSampleSequence()
{
var n = new InverseGaussian(1.0, 2.0);
var ied = n.Samples();
GC.KeepAlive(ied.Take(5).ToArray());
}
/// <summary>
/// Can estimate distribution parameters.
/// </summary>
[TestCase(10.0, 0.1)]
[TestCase(1.5, 1.5)]
[TestCase(2.5, 2.5)]
[TestCase(2.5, 5.0)]
[TestCase(10.0, 50.0)]
public void CanEstimateParameters(double mu, double lambda)
{
var original = new InverseGaussian(mu, lambda, new Numerics.Random.MersenneTwister(100));
var estimated = InverseGaussian.Estimate(original.Samples().Take(1000000));
AssertHelpers.AlmostEqualRelative(mu, estimated.Mu, 1);
AssertHelpers.AlmostEqualRelative(lambda, estimated.Lambda, 1);
}
}
}

446
src/Numerics.Tests/DistributionTests/Continuous/TruncatedParetoTests.cs

@ -0,0 +1,446 @@
// <copyright file="TruncatedParetoTests.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
//
// Copyright (c) 2009-2019 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>
using MathNet.Numerics.Distributions;
using NUnit.Framework;
using System;
using System.Linq;
namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
{
/// <summary>
/// <c>TruncatedPareto</c> distribution tests.
/// </summary>
[TestFixture, Category("Distributions")]
public class TruncatedParetoTests
{
private const int highPrecision = 12;
private const int lowPrecision = 6;
private const double highTruncation = 1E8;
/// <summary>
/// Can create <c>TruncatedPareto</c>.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// /// <param name="truncation">Truncation parameter.</param>
[TestCase(2.0, 2.0, 500.0)]
[TestCase(10.0, 2.0, 1000.0)]
[TestCase(100.0, 4.0, 10000.0)]
[TestCase(3.0, 10.0, 5.0)]
[TestCase(10.0, 100.0, 20.0)]
public void CanCreateTruncatedPareto(double scale, double shape, double truncation)
{
var n = new TruncatedPareto(scale, shape, truncation);
Assert.AreEqual(scale, n.Scale, highPrecision);
Assert.AreEqual(shape, n.Shape, highPrecision);
Assert.AreEqual(truncation, n.Truncation, highPrecision);
}
/// <summary>
/// Can create <c>TruncatedPareto</c> with random source.
/// </summary>
[Test]
public void CanCreateTruncatedParetoWithRandomSource()
{
var randomSource = new Numerics.Random.MersenneTwister(100);
var n = new TruncatedPareto(10.0, 10.0, 1000, randomSource);
Assert.AreEqual(randomSource, n.RandomSource);
}
/// <summary>
/// <c>TruncatedPareto</c> create fails with bad parameters.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// /// <param name="truncation">Truncation parameter.</param>
[TestCase(1.0, -1.0, 15.0)]
[TestCase(1.0, 3.0, 0.5)]
[TestCase(-1.0, 2.0, 15.0)]
[TestCase(double.NaN, 1.0, 1.0)]
[TestCase(1.0, double.PositiveInfinity, 0.0)]
[TestCase(5.0, 2.0, double.PositiveInfinity)]
public void TruncatedParetoCreateFailsWithBadParameters(double scale, double shape, double truncation)
{
Assert.That(() => new TruncatedPareto(scale, shape, truncation), Throws.ArgumentException);
}
/// <summary>
/// Validate IsValidParameterSet.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// /// <param name="truncation">Truncation parameter.</param>
[TestCase(2.0, 2.0, 500, true)]
[TestCase(10.0, 10.0, 1000, true)]
[TestCase(5.0, 3.0, 1, false)]
public void ValidateIsValidParameterSet(double scale, double shape, double truncation, bool validity)
{
Assert.AreEqual(TruncatedPareto.IsValidParameterSet(scale, shape, truncation), validity);
}
/// <summary>
/// Validate ToString.
/// </summary>
[Test]
public void ValidateToString()
{
var n = new TruncatedPareto(1d, 2d, 100d);
Assert.AreEqual("Truncated Pareto(Scale = 1, Shape = 2, Truncation = 100)", n.ToString());
}
/// <summary>
/// Validate mode throws not supported exception
/// </summary>
/// <param name="scale">Scale value.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
[TestCase(10.0, 10.0, highTruncation)]
public void ValidateModeFailsWithNotSupported(double scale, double shape, double truncation)
{
var n = new TruncatedPareto(scale, shape, truncation);
Assert.Throws<NotSupportedException>(() => { var x = n.Mode; });
}
/// <summary>
/// Validate median.
/// </summary>
/// <param name="scale">Scale value.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="median">Expected value.</param>
[TestCase(10.0, 10.0, highTruncation, 10.717734625362931)]
[TestCase(100, 3.5, 10000, 1.219013619375516e+02)]
[TestCase(1000, 5.5, 100000, 1.134312522193400e+03)]
public void ValidateMedian(double scale, double shape, double truncation, double median)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(median, n.Median, highPrecision);
}
/// <summary>
/// Validate mean.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="std">Expected value.</param>
[TestCase(10.0, 10.0, highTruncation, 11.111111111111111)]
[TestCase(100, 3.5, 10000, 1.399986139998614e+02)]
public void ValidateMean(double scale, double shape, double truncation, double mean)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(mean, n.Mean, highPrecision);
}
/// <summary>
/// Validate variance.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="variance">Expected value.</param>
[TestCase(10.0, 10.0, highTruncation, 1.543209876543210)]
[TestCase(100, 3.5, 10000, 3.710390409118045e+03)]
[TestCase(1000, 5.5, 100000, 7.760125676537910e+04)]
public void ValidateVariance(double scale, double shape, double truncation, double variance)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(variance, n.Variance, highPrecision);
}
/// <summary>
/// Validate standard deviation.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="std">Expected value.</param>
[TestCase(10.0, 10.0, highTruncation, 1.242259987499883)]
[TestCase(100, 3.5, 10000, 60.912974062329653)]
[TestCase(1000, 5.5, 100000, 2.785700212969427e+02)]
public void ValidateStandardDeviation(double scale, double shape, double truncation, double std)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(std, n.StdDev, highPrecision);
}
/// <summary>
/// Validate minimum.
/// </summary>
[Test]
public void ValidateMinimum()
{
var n = new TruncatedPareto(1.0, 2.0, 500.0);
Assert.AreEqual(1.0, n.Minimum);
}
/// <summary>
/// Validate maximum.
/// </summary>
[Test]
public void ValidateMaximum()
{
var n = new TruncatedPareto(1.0, 2.0, 500.0);
Assert.AreEqual(500.0, n.Maximum);
}
/// <summary>
/// Validate entropy throws not supported exception
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="lambda">Lambda parameter.</param>
[TestCase(100, 3.5, 10000)]
public void ValidateEntropyFailsWithNotSupported(double scale, double shape, double truncation)
{
var n = new TruncatedPareto(scale, shape, truncation);
Assert.Throws<NotSupportedException>(() => { var x = n.Entropy; });
}
/// <summary>
/// Validate GetMoments.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="n"> order of the moment.</param>
/// <param name="value">Expected value.</param>
[TestCase(10.0, 10.0, highTruncation, 1, 11.111111111111111)]
[TestCase(1.0, 6.0, 10000, 1, 1.2)]
[TestCase(100, 3.5, 10000, 1, 1.399986139998614e+02)]
[TestCase(100, 3.5, 10000, 2, 2.331000233100023e+04)]
[TestCase(100, 3.5, 10000, 3, 6.300000630000063e+06)]
[TestCase(1000, 5.5, 100000, 1, 1.222222221012222e+03)]
[TestCase(1000, 5.5, 100000, 2, 1.571428414301428e+06)]
[TestCase(1000, 5.5, 100000, 4, 3.663000000036630e+12)]
[TestCase(1000, 2.0, 100000, 2, 9.211261498125996e+06)]
[TestCase(1000, 3.0, 100000, 3, 1.381552437348865e+10)]
public void ValidateMoments(double scale, double shape, double truncation, int order, double value)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(value, n.GetMoment(order), highPrecision);
}
/// <summary>
/// Validate skewness.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="skewness">Expected value.</param>
[TestCase(10.0, 10.0, highTruncation, 2.8110568859997356)]
public void ValidateSkewness(double scale, double shape, double truncation, double skewness)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(skewness, n.Skewness, highPrecision);
}
/// <summary>
/// Validate density.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="x">The location at which to compute the density.</param>
/// <param name="p">Expected value.</param>
[TestCase(1, 1, highTruncation, 1, 1)]
[TestCase(1, 1, highTruncation, 1.5, 4 / 9.0)]
[TestCase(1, 1, highTruncation, 5, 1 / 25.0)]
[TestCase(1, 1, highTruncation, 50, 1 / 2500.0)]
[TestCase(1, 4, highTruncation, 1, 4)]
[TestCase(1, 4, highTruncation, 1.5, 128 / 243.0)]
[TestCase(1, 4, highTruncation, 50, 1 / 78125000.0)]
[TestCase(3, 2, highTruncation, 3, 2 / 3.0)]
[TestCase(3, 2, highTruncation, 5, 18 / 125.0)]
[TestCase(25, 100, highTruncation, 50, 1.5777218104420236e-30)]
[TestCase(100, 25, highTruncation, 150, 6.6003546737276816e-6)]
[TestCase(100, 3.5, 10000, 1000, 1.106797291738662e-06)]
[TestCase(100, 3.5, 10000, 2000, 4.891399189920708e-08)]
[TestCase(1000, 5.5, 100000, 2000, 6.076698900882660e-05)]
public void ValidateDensity(double scale, double shape, double truncation, double x, double p)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(p, n.Density(x), lowPrecision);
AssertHelpers.AlmostEqualRelative(p, TruncatedPareto.PDF(scale, shape, truncation, x), lowPrecision);
}
/// <summary>
/// Validate density log.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="x">The location at which to compute the log density.</param>
/// <param name="p">Expected value.</param>
[TestCase(3, 2, highTruncation, 3, -0.405465108108164)]
[TestCase(100, 3.5, 10000, 1000, -13.714040035965924)]
[TestCase(100, 3.5, 10000, 2000, -16.833202348485678)]
public void ValidateDensityLn(double scale, double shape, double truncation, double x, double p)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(p, n.DensityLn(x), lowPrecision);
AssertHelpers.AlmostEqualRelative(p, TruncatedPareto.PDFLn(scale, shape, truncation, x), lowPrecision);
}
/// <summary>
/// Validate cumulative distribution.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="x">The location at which to compute the cumulative distribution function.</param>
/// <param name="expected">Expected value.</param>
[TestCase(1.0, 1.0, highTruncation, 1.0, 0)]
[TestCase(7.0, 7.0, highTruncation, 10.0, 0.9176457)]
[TestCase(10.0, 10.0, highTruncation, 12.0, 0.83849441711015427)]
[TestCase(100, 3.5, 10000, 102, 0.066961862452387)]
[TestCase(100, 3.5, 10000, 1000, 0.999683872202370)]
[TestCase(100, 3.5, 10000, 2000, 0.999972149147496)]
[TestCase(1000, 5.5, 100000, 2000, 0.977902913097699)]
public void ValidateCumulativeDistribution(double scale, double shape, double truncation, double x, double expected)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(expected, n.CumulativeDistribution(x), highPrecision);
AssertHelpers.AlmostEqualRelative(expected, TruncatedPareto.CDF(scale, shape, truncation, x), highPrecision);
}
/// <summary>
/// Validate inverse cumulative distribution.
/// </summary>
/// <param name="scale">Scale parameter.</param>
/// <param name="shape">Shape parameter.</param>
/// <param name="truncation">Truncation parameter.</param>
/// <param name="p">The location at which to compute the inverse cumulative distribution function.</param>
/// <param name="expected">Expected value.</param>
[TestCase(1.0, 1.0, highTruncation, 0, 1.0)]
[TestCase(7.0, 7.0, highTruncation, 0.9176457, 10.0)]
[TestCase(10.0, 10.0, highTruncation, 0.83849441711015427, 12.0)]
[TestCase(100, 3.5, 10000, 0.066961862452387, 102)]
[TestCase(100, 3.5, 10000, 0.999683872202370, 1000)]
[TestCase(100, 3.5, 10000, 0.999972149147496, 2000)]
[TestCase(1000, 5.5, 100000, 0.977902913097699, 2000)]
public void ValidateInverseCumulativeDistribution(double scale, double shape, double truncation, double p, double expected)
{
var n = new TruncatedPareto(scale, shape, truncation);
AssertHelpers.AlmostEqualRelative(expected, n.InvCDF(p), lowPrecision);
AssertHelpers.AlmostEqualRelative(expected, TruncatedPareto.ICDF(scale, shape, truncation, p), lowPrecision);
}
/// <summary>
/// Can sample static.
/// </summary>
[Test]
public void CanSampleStatic()
{
TruncatedPareto.Sample(new Numerics.Random.MersenneTwister(100), 10.0, 10.0, 1000.0);
}
/// <summary>
/// Can fill array with samples static.
/// </summary>
[Test]
public void CanFillSampleArrayStatic()
{
double[] samples = new double[100];
TruncatedPareto.Samples(new Numerics.Random.MersenneTwister(100), samples, 10.0, 10.0, 1000.0);
Assert.IsTrue(!samples.Any(x => x == 0));
}
/// <summary>
/// Can sample sequence static.
/// </summary>
[Test]
public void CanSampleSequenceStatic()
{
var ied = TruncatedPareto.Samples(new Numerics.Random.MersenneTwister(100), 10.0, 10.0, 1000.0);
GC.KeepAlive(ied.Take(5).ToArray());
}
/// <summary>
/// Fail sample static with bad parameters.
/// </summary>
[Test]
public void FailSampleStatic()
{
Assert.That(() => { var d = TruncatedPareto.Sample(new Numerics.Random.MersenneTwister(100), 10.0, 10.0, 5.0); }, Throws.ArgumentException);
}
/// <summary>
/// Fail filling array with samples with bad parameters static.
/// </summary>
[Test]
public void FailFillingSampleArrayStatic()
{
double[] samples = new double[100];
Assert.That(() => { TruncatedPareto.Samples(new Numerics.Random.MersenneTwister(100), samples, 10.0, 10.0, 5.0); }, Throws.ArgumentException);
}
/// <summary>
/// Fail sample sequence static with bad parameters.
/// </summary>
[Test]
public void FailSampleSequenceStatic()
{
Assert.That(() => { var ied = TruncatedPareto.Samples(new Numerics.Random.MersenneTwister(100), 10.0, 10.0, 5.0).First(); }, Throws.ArgumentException);
}
/// <summary>
/// Can sample.
/// </summary>
[Test]
public void CanSample()
{
var n = new TruncatedPareto(10.0, 10.0, 1000.0);
n.Sample();
}
/// <summary>
/// Can fill array with samples.
/// </summary>
[Test]
public void CanFillSampleArray()
{
double[] samples = new double[100];
var n = new TruncatedPareto(10.0, 10.0, 1000.0);
n.Samples(samples);
Assert.IsTrue(!samples.Any(x => x == 0));
}
/// <summary>
/// Can sample sequence.
/// </summary>
[Test]
public void CanSampleSequence()
{
var n = new TruncatedPareto(10.0, 10.0, 1000.0);
var ied = n.Samples();
GC.KeepAlive(ied.Take(5).ToArray());
}
}
}

429
src/Numerics/Distributions/Burr.cs

@ -0,0 +1,429 @@
// <copyright file="Burr.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
//
// Copyright (c) 2009-2019 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>
using MathNet.Numerics.Properties;
using MathNet.Numerics.Random;
using System;
using System.Collections.Generic;
namespace MathNet.Numerics.Distributions
{
public class Burr : IContinuousDistribution
{
private System.Random _random;
/// <summary>
/// Gets the scale (a) of the distribution. Range: a > 0.
/// </summary>
public double a { get; }
/// <summary>
/// Gets the first shape parameter (c) of the distribution. Range: c > 0.
/// </summary>
public double c { get; }
/// <summary>
/// Gets the second shape parameter (k) of the distribution. Range: k > 0.
/// </summary>
public double k { get; }
/// <summary>
/// Initializes a new instance of the Burr Type XII class.
/// </summary>
/// <param name="a">The scale parameter a of the Burr distribution. Range: a > 0.</param>
/// <param name="c">The first shape parameter c of the Burr distribution. Range: c > 0.</param>
/// <param name="k">The second shape parameter k of the Burr distribution. Range: k > 0.</param>
/// <param name="randomSource">The random number generator which is used to draw random samples. Optional, can be null.</param>
public Burr(double a, double c, double k, System.Random randomSource = null)
{
if (!IsValidParameterSet(a, c, k))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
_random = randomSource ?? SystemRandomSource.Default;
this.a = a;
this.c = c;
this.k = k;
}
/// <summary>
/// A string representation of the distribution.
/// </summary>
/// <returns>a string representation of the distribution.</returns>
public override string ToString()
{
return "Burr(a = " + a + ", c = " + c + ", k = " + k + ")";
}
/// <summary>
/// Tests whether the provided values are valid parameters for this distribution.
/// </summary>
/// <param name="a">The scale parameter a of the Burr distribution. Range: a > 0.</param>
/// <param name="c">The first shape parameter c of the Burr distribution. Range: c > 0.</param>
/// <param name="k">The second shape parameter k of the Burr distribution. Range: k > 0.</param>
public static bool IsValidParameterSet(double a, double c, double k)
{
var allFinite = a.IsFinite() && c.IsFinite() && k.IsFinite();
return allFinite && a > 0.0 && c > 0.0 && k > 0.0;
}
/// <summary>
/// Gets the random number generator which is used to draw random samples.
/// </summary>
public System.Random RandomSource
{
get { return _random; }
set { _random = value ?? SystemRandomSource.Default; }
}
/// <summary>
/// Gets the mean of the Burr distribution.
/// </summary>
public double Mean
{
get
{
return (1 / SpecialFunctions.Gamma(k)) * a * SpecialFunctions.Gamma(1 + 1 / c) * SpecialFunctions.Gamma(k - 1 / c);
}
}
/// <summary>
/// Gets the variance of the Burr distribution.
/// </summary>
public double Variance
{
get
{
return (1 / SpecialFunctions.Gamma(k)) * Math.Pow(a, 2) * SpecialFunctions.Gamma(1 + 2 / c) * SpecialFunctions.Gamma(k - 2 / c)
- Math.Pow((1 / SpecialFunctions.Gamma(k)) * a * SpecialFunctions.Gamma(1 + 1 / c) * SpecialFunctions.Gamma(k - 1 / c), 2);
}
}
/// <summary>
/// Gets the standard deviation of the Burr distribution.
/// </summary>
public double StdDev
{
get
{
return Math.Sqrt(Variance);
}
}
/// <summary>
/// Gets the mode of the Burr distribution.
/// </summary>
public double Mode
{
get
{
return a * Math.Pow((c - 1) / (c * k + 1), 1 / c);
}
}
/// <summary>
/// Gets the minimum of the Burr distribution.
/// </summary>
public double Minimum
{
get { return 0.0; }
}
/// <summary>
/// Gets the maximum of the Burr distribution.
/// </summary>
public double Maximum
{
get { return double.PositiveInfinity; }
}
/// <summary>
/// Gets the entropy of the Burr distribution (currently not supported).
/// </summary>
public double Entropy
{
get { throw new NotSupportedException(); }
}
/// <summary>
/// Gets the skewness of the Burr distribution.
/// </summary>
public double Skewness
{
get
{
var mean = Mean;
var variance = Variance;
var std = StdDev;
return (GetMoment(3) - 3 * mean * variance - mean * mean * mean) / (std * std * std);
}
}
/// <summary>
/// Gets the median of the Burr distribution.
/// </summary>
public double Median
{
get
{
return a * Math.Pow(Math.Pow(2, 1 / k) - 1, 1 / c);
}
}
/// <summary>
/// Generates a sample from the Burr distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return SampleUnchecked(_random, a, c, k);
}
/// <summary>
/// Fills an array with samples generated from the distribution.
/// </summary>
/// <param name="values">The array to fill with the samples.</param>
public void Samples(double[] values)
{
SamplesUnchecked(_random, values, a, c, k);
}
/// <summary>
/// Generates a sequence of samples from the Burr distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
return SamplesUnchecked(_random, a, c, k);
}
/// <summary>
/// Generates a sample from the Burr distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="a">The scale parameter a of the Burr distribution. Range: a > 0.</param>
/// <param name="c">The first shape parameter c of the Burr distribution. Range: c > 0.</param>
/// <param name="k">The second shape parameter k of the Burr distribution. Range: k > 0.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(System.Random rnd, double a, double c, double k)
{
if (!IsValidParameterSet(a, c, k))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, a, c, k);
}
/// <summary>
/// Fills an array with samples generated from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="values">The array to fill with the samples.</param>
/// <param name="a">The scale parameter a of the Burr distribution. Range: a > 0.</param>
/// <param name="c">The first shape parameter c of the Burr distribution. Range: c > 0.</param>
/// <param name="k">The second shape parameter k of the Burr distribution. Range: k > 0.</param>
public static void Samples(System.Random rnd, double[] values, double a, double c, double k)
{
if (!IsValidParameterSet(a, c, k))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
SamplesUnchecked(rnd, values, a, c, k);
}
/// <summary>
/// Generates a sequence of samples from the Burr distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="a">The scale parameter a of the Burr distribution. Range: a > 0.</param>
/// <param name="c">The first shape parameter c of the Burr distribution. Range: c > 0.</param>
/// <param name="k">The second shape parameter k of the Burr distribution. Range: k > 0.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(System.Random rnd, double a, double c, double k)
{
if (!IsValidParameterSet(a, c, k))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return SamplesUnchecked(rnd, a, c, k);
}
internal static double SampleUnchecked(System.Random rnd, double a, double c, double k)
{
var k_inv = 1 / k;
var c_inv = 1 / c;
double u = rnd.NextDouble();
return a * Math.Pow(Math.Pow(1 - u, -k_inv) - 1, c_inv);
}
internal static void SamplesUnchecked(System.Random rnd, double[] values, double a, double c, double k)
{
if (values.Length == 0)
{
return;
}
var k_inv = 1 / k;
var c_inv = 1 / c;
double[] u = rnd.NextDoubles(values.Length);
for (var j = 0; j < values.Length; ++j)
{
values[j] = a * Math.Pow(Math.Pow(1 - u[j], -k_inv) - 1, c_inv);
}
}
internal static IEnumerable<double> SamplesUnchecked(System.Random rnd, double a, double c, double k)
{
while (true)
{
yield return SampleUnchecked(rnd, a, c, k);
}
}
/// <summary>
/// Gets the n-th raw moment of the distribution.
/// </summary>
/// <param name="n">The order (n) of the moment. Range: n ≥ 1.</param>
/// <returns>the n-th moment of the distribution.</returns>
public double GetMoment(double n)
{
if (n > k * c)
{
throw new ArgumentException(Resources.ArgumentParameterSetInvalid);
}
var lambda_n = (n / c) * SpecialFunctions.Gamma(n / c) * SpecialFunctions.Gamma(k - n / c);
return Math.Pow(a, n) * lambda_n / SpecialFunctions.Gamma(k);
}
/// <summary>
/// Computes the probability density of the distribution (PDF) at x, i.e. ∂P(X ≤ x)/∂x.
/// </summary>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
/// <seealso cref="PDF"/>
public double Density(double x)
{
return DensityImpl(a, c, k, x);
}
/// <summary>
/// Computes the log probability density of the distribution (lnPDF) at x, i.e. ln(∂P(X ≤ x)/∂x).
/// </summary>
/// <param name="x">The location at which to compute the log density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
/// <seealso cref="PDFLn"/>
public double DensityLn(double x)
{
return DensityLnImpl(a, c, k, x);
}
/// <summary>
/// Computes the cumulative distribution (CDF) of the distribution at x, i.e. P(X ≤ x).
/// </summary>
/// <param name="x">The location at which to compute the cumulative distribution function.</param>
/// <returns>the cumulative distribution at location <paramref name="x"/>.</returns>
/// <seealso cref="CDF"/>
public double CumulativeDistribution(double x)
{
return CumulativeDistributionImpl(a, c, k, x);
}
/// <summary>
/// Computes the probability density of the distribution (PDF) at x, i.e. ∂P(X ≤ x)/∂x.
/// </summary>
/// <param name="a">The scale parameter a of the Burr distribution. Range: a > 0.</param>
/// <param name="c">The first shape parameter c of the Burr distribution. Range: c > 0.</param>
/// <param name="k">The second shape parameter k of the Burr distribution. Range: k > 0.</param>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
/// <seealso cref="Density"/>
public static double PDF(double a, double c, double k, double x)
{
if (!IsValidParameterSet(a, c, k))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return DensityImpl(a, c, k, x);
}
/// <summary>
/// Computes the log probability density of the distribution (lnPDF) at x, i.e. ln(∂P(X ≤ x)/∂x).
/// </summary>
/// <param name="a">The scale parameter a of the Burr distribution. Range: a > 0.</param>
/// <param name="c">The first shape parameter c of the Burr distribution. Range: c > 0.</param>
/// <param name="k">The second shape parameter k of the Burr distribution. Range: k > 0.</param>
/// <param name="x">The location at which to compute the log density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
/// <seealso cref="DensityLn"/>
public static double PDFLn(double a, double c, double k, double x)
{
if (!IsValidParameterSet(a, c, k))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return DensityLnImpl(a, c, k, x);
}
/// <summary>
/// Computes the cumulative distribution (CDF) of the distribution at x, i.e. P(X ≤ x).
/// </summary>
/// <param name="a">The scale parameter a of the Burr distribution. Range: a > 0.</param>
/// <param name="c">The first shape parameter c of the Burr distribution. Range: c > 0.</param>
/// <param name="k">The second shape parameter k of the Burr distribution. Range: k > 0.</param>
/// <param name="x">The location at which to compute the cumulative distribution function.</param>
/// <returns>the cumulative distribution at location <paramref name="x"/>.</returns>
/// <seealso cref="CumulativeDistribution"/>
public static double CDF(double a, double c, double k, double x)
{
if (!IsValidParameterSet(a, c, k))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return CumulativeDistributionImpl(a, c, k, x);
}
internal static double DensityImpl(double a, double c, double k, double x)
{
var numerator = (k * c / a) * Math.Pow(x / a, c - 1);
var denominator = Math.Pow(1 + Math.Pow(x / a, c), k + 1);
return numerator / denominator;
}
internal static double DensityLnImpl(double a, double c, double k, double x)
{
return Math.Log(DensityImpl(a, c, k, x));
}
internal static double CumulativeDistributionImpl(double a, double c, double k, double x)
{
var denominator = Math.Pow(1 + Math.Pow(x / a, c), k);
return 1 - 1 / denominator;
}
}
}

449
src/Numerics/Distributions/InverseGaussian.cs

@ -0,0 +1,449 @@
// <copyright file="InverseGaussian.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
//
// Copyright (c) 2009-2019 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>
using MathNet.Numerics.Properties;
using MathNet.Numerics.Random;
using MathNet.Numerics.Statistics;
using System;
using System.Collections.Generic;
using System.Linq;
namespace MathNet.Numerics.Distributions
{
public class InverseGaussian : IContinuousDistribution
{
private System.Random _random;
/// <summary>
/// Gets the mean (μ) of the distribution. Range: μ > 0.
/// </summary>
public double Mu { get; }
/// <summary>
/// Gets the shape (λ) of the distribution. Range: λ > 0.
/// </summary>
public double Lambda { get; }
/// <summary>
/// Initializes a new instance of the InverseGaussian class.
/// </summary>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
/// <param name="randomSource">The random number generator which is used to draw random samples.</param>
public InverseGaussian(double mu, double lambda, System.Random randomSource = null)
{
if (!IsValidParameterSet(mu, lambda))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
_random = randomSource ?? SystemRandomSource.Default;
Mu = mu;
Lambda = lambda;
}
/// <summary>
/// A string representation of the distribution.
/// </summary>
/// <returns>a string representation of the distribution.</returns>
public override string ToString()
{
return "InverseGaussian(μ = " + Mu + ", λ = " + Lambda + ")";
}
/// <summary>
/// Tests whether the provided values are valid parameters for this distribution.
/// </summary>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
public static bool IsValidParameterSet(double mu, double lambda)
{
var allFinite = mu.IsFinite() && lambda.IsFinite();
return allFinite && mu > 0.0 && lambda > 0.0;
}
/// <summary>
/// Gets the random number generator which is used to draw random samples.
/// </summary>
public System.Random RandomSource
{
get { return _random; }
set { _random = value ?? SystemRandomSource.Default; }
}
/// <summary>
/// Gets the mean of the Inverse Gaussian distribution.
/// </summary>
public double Mean
{
get
{
return Mu;
}
}
/// <summary>
/// Gets the variance of the Inverse Gaussian distribution.
/// </summary>
public double Variance
{
get
{
return Math.Pow(Mu, 3) / Lambda;
}
}
/// <summary>
/// Gets the standard deviation of the Inverse Gaussian distribution.
/// </summary>
public double StdDev
{
get
{
return Math.Sqrt(Variance);
}
}
/// <summary>
/// Gets the median of the Inverse Gaussian distribution.
/// No closed form analytical expression exists, so this value is approximated numerically and can throw an exception.
/// </summary>
public double Median
{
get { return InvCDF(0.5); }
}
/// <summary>
/// Gets the minimum of the Inverse Gaussian distribution.
/// </summary>
public double Minimum
{
get { return 0.0; }
}
/// <summary>
/// Gets the maximum of the Inverse Gaussian distribution.
/// </summary>
public double Maximum
{
get { return double.PositiveInfinity; }
}
/// <summary>
/// Gets the skewness of the Inverse Gaussian distribution.
/// </summary>
public double Skewness
{
get { return 3 * Math.Sqrt(Mu / Lambda); }
}
/// <summary>
/// Gets the kurtosis of the Inverse Gaussian distribution.
/// </summary>
public double Kurtosis
{
get { return 15 * Mu / Lambda; }
}
/// <summary>
/// Gets the mode of the Inverse Gaussian distribution.
/// </summary>
public double Mode
{
get { return Mu * (Math.Sqrt(1 + (9 * Mu * Mu) / (4 * Lambda * Lambda)) - (3 * Mu) / (2 * Lambda)); }
}
/// <summary>
/// Gets the entropy of the Inverse Gaussian distribution (currently not supported).
/// </summary>
public double Entropy
{
get { throw new NotSupportedException(); }
}
/// <summary>
/// Generates a sample from the inverse Gaussian distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return SampleUnchecked(_random, Mu, Lambda);
}
/// <summary>
/// Fills an array with samples generated from the distribution.
/// </summary>
/// <param name="values">The array to fill with the samples.</param>
public void Samples(double[] values)
{
SamplesUnchecked(_random, values, Mu, Lambda);
}
/// <summary>
/// Generates a sequence of samples from the inverse Gaussian distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
return SamplesUnchecked(_random, Mu, Lambda);
}
/// <summary>
/// Generates a sample from the inverse Gaussian distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(System.Random rnd, double mu, double lambda)
{
if (!IsValidParameterSet(mu, lambda))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, mu, lambda);
}
/// <summary>
/// Fills an array with samples generated from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="values">The array to fill with the samples.</param>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
public static void Samples(System.Random rnd, double[] values, double mu, double lambda)
{
if (!IsValidParameterSet(mu, lambda))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
SamplesUnchecked(rnd, values, mu, lambda);
}
/// <summary>
/// Generates a sequence of samples from the Burr distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(System.Random rnd, double mu, double lambda)
{
if (!IsValidParameterSet(mu, lambda))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return SamplesUnchecked(rnd, mu, lambda);
}
internal static double SampleUnchecked(System.Random rnd, double mu, double lambda)
{
double v = MathNet.Numerics.Distributions.Normal.Sample(rnd, 0, 1);
double test = rnd.NextDouble();
return InverseGaussianSampleImpl(mu, lambda, v, test);
}
internal static void SamplesUnchecked(System.Random rnd, double[] values, double mu, double lambda)
{
if (values.Length == 0)
{
return;
}
double[] v = new double[values.Length];
MathNet.Numerics.Distributions.Normal.Samples(rnd, v, 0, 1);
double[] test = rnd.NextDoubles(values.Length);
for (var j = 0; j < values.Length; ++j)
{
values[j] = InverseGaussianSampleImpl(mu, lambda, v[j], test[j]);
}
}
internal static IEnumerable<double> SamplesUnchecked(System.Random rnd, double mu, double lambda)
{
while (true)
{
yield return SampleUnchecked(rnd, mu, lambda);
}
}
internal static double InverseGaussianSampleImpl(double mu, double lambda, double normalSample, double uniformSample)
{
double y = normalSample * normalSample;
double x = mu + (mu * mu * y) / (2 * lambda) - (mu / (2 * lambda)) * Math.Sqrt(4 * mu * lambda * y + mu * mu * y * y);
if (uniformSample <= mu / (mu + x))
return x;
else
return mu * mu / x;
}
/// <summary>
/// Computes the probability density of the distribution (PDF) at x, i.e. ∂P(X ≤ x)/∂x.
/// </summary>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
/// <seealso cref="PDF"/>
public double Density(double x)
{
return DensityImpl(Mu, Lambda, x);
}
/// <summary>
/// Computes the log probability density of the distribution (lnPDF) at x, i.e. ln(∂P(X ≤ x)/∂x).
/// </summary>
/// <param name="x">The location at which to compute the log density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
/// <seealso cref="PDFLn"/>
public double DensityLn(double x)
{
return DensityLnImpl(Mu, Lambda, x);
}
/// <summary>
/// Computes the cumulative distribution (CDF) of the distribution at x, i.e. P(X ≤ x).
/// </summary>
/// <param name="x">The location at which to compute the cumulative distribution function.</param>
/// <returns>the cumulative distribution at location <paramref name="x"/>.</returns>
/// <seealso cref="CDF"/>
public double CumulativeDistribution(double x)
{
return CumulativeDistributionImpl(Mu, Lambda, x);
}
/// <summary>
/// Computes the inverse cumulative distribution (CDF) of the distribution at p, i.e. solving for P(X ≤ x) = p.
/// </summary>
/// <param name="p">The location at which to compute the inverse cumulative distribution function.</param>
/// <returns>the inverse cumulative distribution at location <paramref name="p"/>.</returns>
public double InvCDF(double p)
{
Func<double, double> equationToSolve = (x) => CumulativeDistribution(x) - p;
if (RootFinding.NewtonRaphson.TryFindRoot(equationToSolve, Density, Mode, 0, double.PositiveInfinity, 1e-8, 100, out double quantile))
return quantile;
else
throw new NonConvergenceException(Resources.NumericalEstimationFailed);
}
/// <summary>
/// Computes the probability density of the distribution (PDF) at x, i.e. ∂P(X ≤ x)/∂x.
/// </summary>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
/// <seealso cref="Density"/>
public static double PDF(double mu, double lambda, double x)
{
if (!IsValidParameterSet(mu, lambda))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return DensityImpl(mu, lambda, x);
}
/// <summary>
/// Computes the log probability density of the distribution (lnPDF) at x, i.e. ln(∂P(X ≤ x)/∂x).
/// </summary>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
/// <param name="x">The location at which to compute the log density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
/// <seealso cref="DensityLn"/>
public static double PDFLn(double mu, double lambda, double x)
{
if (!IsValidParameterSet(mu, lambda))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return DensityLnImpl(mu, lambda, x);
}
/// <summary>
/// Computes the cumulative distribution (CDF) of the distribution at x, i.e. P(X ≤ x).
/// </summary>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
/// <param name="x">The location at which to compute the cumulative distribution function.</param>
/// <returns>the cumulative distribution at location <paramref name="x"/>.</returns>
/// <seealso cref="CumulativeDistribution"/>
public static double CDF(double mu, double lambda, double x)
{
if (!IsValidParameterSet(mu, lambda))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return CumulativeDistributionImpl(mu, lambda, x);
}
/// <summary>
/// Computes the inverse cumulative distribution (CDF) of the distribution at p, i.e. solving for P(X ≤ x) = p.
/// </summary>
/// <param name="mu">The mean (μ) of the distribution. Range: μ > 0.</param>
/// <param name="lambda">The shape (λ) of the distribution. Range: λ > 0.</param>
/// <param name="p">The location at which to compute the inverse cumulative distribution function.</param>
/// <returns>the inverse cumulative distribution at location <paramref name="p"/>.</returns>
/// <seealso cref="CumulativeDistribution"/>
public static double ICDF(double mu, double lambda, double p)
{
if (!IsValidParameterSet(mu, lambda))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
var igDist = new InverseGaussian(mu, lambda);
return igDist.InvCDF(p);
}
/// <summary>
/// Estimates the Inverse Gaussian parameters from sample data with maximum-likelihood.
/// </summary>
/// <param name="samples">The samples to estimate the distribution parameters from.</param>
/// <param name="randomSource">The random number generator which is used to draw random samples. Optional, can be null.</param>
/// <returns>An Inverse Gaussian distribution.</returns>
public static InverseGaussian Estimate(IEnumerable<double> samples, System.Random randomSource = null)
{
var sampleVec = samples.ToArray();
var muHat = sampleVec.Mean();
var lambdahat = 1 / (1 / samples.HarmonicMean() - 1 / muHat);
return new InverseGaussian(muHat, lambdahat, randomSource);
}
internal static double DensityImpl(double mu, double lambda, double x)
{
return Math.Sqrt(lambda / (2 * Math.PI * Math.Pow(x, 3))) * Math.Exp(-((lambda * Math.Pow(x - mu, 2)) / (2 * mu * mu * x)));
}
internal static double DensityLnImpl(double mu, double lambda, double x)
{
return Math.Log(Math.Sqrt(lambda / (2 * Math.PI * Math.Pow(x, 3)))) - ((lambda * Math.Pow(x - mu, 2)) / (2 * mu * mu * x));
}
internal static double CumulativeDistributionImpl(double mu, double lambda, double x)
{
return Normal.CDF(0, 1, Math.Sqrt(lambda / x) * (x / mu - 1)) + Math.Exp(2 * lambda / mu) * Normal.CDF(0, 1, -Math.Sqrt(lambda / x) * (x / mu + 1));
}
}
}

466
src/Numerics/Distributions/TruncatedPareto.cs

@ -0,0 +1,466 @@
// <copyright file="TruncatedPareto.cs" company="Math.NET">
// Math.NET Numerics, part of the Math.NET Project
// http://numerics.mathdotnet.com
// http://github.com/mathnet/mathnet-numerics
//
// Copyright (c) 2009-2019 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>
using MathNet.Numerics.Properties;
using MathNet.Numerics.Random;
using System;
using System.Collections.Generic;
namespace MathNet.Numerics.Distributions
{
public class TruncatedPareto : IContinuousDistribution
{
private System.Random _random;
/// <summary>
/// Initializes a new instance of the TruncatedPareto class.
/// </summary>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
/// <param name="randomSource">The random number generator which is used to draw random samples.</param>
/// <exception cref="ArgumentException">If <paramref name="scale"/> or <paramref name="shape"/> are non-positive or if T ≤ xm.</exception>
public TruncatedPareto(double scale, double shape, double truncation, System.Random randomSource = null)
{
if (!IsValidParameterSet(scale, shape, truncation))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
_random = randomSource ?? SystemRandomSource.Default;
Scale = scale;
Shape = shape;
Truncation = truncation;
}
/// <summary>
/// A string representation of the distribution.
/// </summary>
/// <returns>a string representation of the distribution.</returns>
public override string ToString()
{
return "Truncated Pareto(Scale = " + Scale + ", Shape = " + Shape + ", Truncation = " + Truncation + ")";
}
/// <summary>
/// Tests whether the provided values are valid parameters for this distribution.
/// </summary>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
public static bool IsValidParameterSet(double scale, double shape, double truncation)
{
var allFinite = scale.IsFinite() && shape.IsFinite() && truncation.IsFinite();
return allFinite && scale > 0.0 && shape > 0.0 && truncation > scale;
}
/// <summary>
/// Gets the random number generator which is used to draw random samples.
/// </summary>
public System.Random RandomSource
{
get { return _random; }
set { _random = value ?? SystemRandomSource.Default; }
}
/// <summary>
/// Gets the scale (xm) of the distribution. Range: xm > 0.
/// </summary>
public double Scale { get; }
/// <summary>
/// Gets the shape (α) of the distribution. Range: α > 0.
/// </summary>
public double Shape { get; }
/// <summary>
/// Gets the truncation (T) of the distribution. Range: T > 0.
/// </summary>
public double Truncation { get; }
/// <summary>
/// Gets the n-th raw moment of the distribution.
/// </summary>
/// <param name="n">The order (n) of the moment. Range: n ≥ 1.</param>
/// <returns>the n-th moment of the distribution.</returns>
public double GetMoment(int n)
{
double moment;
if (Shape.AlmostEqual(n))
{
moment = ((Shape * Math.Pow(Scale, n)) / (1 - Math.Pow(Scale / Truncation, Shape))) * Math.Log(Truncation / Scale);
}
else
{
moment = ((Shape * Math.Pow(Scale, n)) / (Shape - n)) * ((1 - Math.Pow((Scale / Truncation), (Shape - n))) / (1 - Math.Pow(Scale / Truncation, Shape)));
}
return moment;
}
/// <summary>
/// Gets the mean of the truncated Pareto distribution.
/// </summary>
public double Mean
{
get
{
return GetMoment(1);
}
}
/// <summary>
/// Gets the variance of the truncated Pareto distribution.
/// </summary>
public double Variance
{
get
{
return GetMoment(2) - Math.Pow(GetMoment(1), 2);
}
}
/// <summary>
/// Gets the standard deviation of the truncated Pareto distribution.
/// </summary>
public double StdDev
{
get
{
return Math.Sqrt(Variance);
}
}
/// <summary>
/// Gets the mode of the truncated Pareto distribution (not supported).
/// </summary>
public double Mode
{
get { throw new NotSupportedException(); }
}
/// <summary>
/// Gets the minimum of the truncated Pareto distribution.
/// </summary>
public double Minimum
{
get { return Scale; }
}
/// <summary>
/// Gets the maximum of the truncated Pareto distribution.
/// </summary>
public double Maximum
{
get { return Truncation; }
}
/// <summary>
/// Gets the entropy of the truncated Pareto distribution (not supported).
/// </summary>
public double Entropy
{
get { throw new NotSupportedException(); }
}
/// <summary>
/// Gets the skewness of the truncated Pareto distribution.
/// </summary>
public double Skewness
{
get
{
var mean = Mean;
var variance = Variance;
var std = StdDev;
return (GetMoment(3) - 3.0 * mean * variance - mean * mean * mean) / (std * std * std);
}
}
/// <summary>
/// Gets the median of the truncated Pareto distribution.
/// </summary>
public double Median
{
get
{
return Scale * Math.Pow(1.0 - (1.0 / 2.0) * (1.0 - Math.Pow(Scale / Truncation, Shape)), -(1.0 / Shape));
}
}
/// <summary>
/// Generates a sample from the truncated Pareto distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return SampleUnchecked(_random, Scale, Shape, Truncation);
}
/// <summary>
/// Fills an array with samples generated from the distribution.
/// </summary>
/// <param name="values">The array to fill with the samples.</param>
public void Samples(double[] values)
{
SamplesUnchecked(_random, values, Scale, Shape, Truncation);
}
/// <summary>
/// Generates a sequence of samples from the truncated Pareto distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
return SamplesUnchecked(_random, Scale, Shape, Truncation);
}
/// <summary>
/// Generates a sample from the truncated Pareto distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(System.Random rnd, double scale, double shape, double truncation)
{
if (!IsValidParameterSet(scale, shape, truncation))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, scale, shape, truncation);
}
/// <summary>
/// Fills an array with samples generated from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="values">The array to fill with the samples.</param>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
public static void Samples(System.Random rnd, double[] values, double scale, double shape, double truncation)
{
if (!IsValidParameterSet(scale, shape, truncation))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
SamplesUnchecked(rnd, values, scale, shape, truncation);
}
/// <summary>
/// Generates a sequence of samples from the truncated Pareto distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(System.Random rnd, double scale, double shape, double truncation)
{
if (!IsValidParameterSet(scale, shape, truncation))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return SamplesUnchecked(rnd, scale, shape, truncation);
}
internal static double SampleUnchecked(System.Random rnd, double scale, double shape, double truncation)
{
double uniform = rnd.NextDouble();
return InvCDFUncheckedImpl(scale, shape, truncation, uniform);
}
internal static void SamplesUnchecked(System.Random rnd, double[] values, double scale, double shape, double truncation)
{
if (values.Length == 0)
{
return;
}
double[] uniforms = rnd.NextDoubles(values.Length);
for (var j = 0; j < values.Length; ++j)
{
values[j] = InvCDFUncheckedImpl(scale, shape, truncation, uniforms[j]);
}
}
internal static IEnumerable<double> SamplesUnchecked(System.Random rnd, double scale, double shape, double truncation)
{
while (true)
{
yield return SampleUnchecked(rnd, scale, shape, truncation);
}
}
/// <summary>
/// Computes the probability density of the distribution (PDF) at x, i.e. ∂P(X ≤ x)/∂x.
/// </summary>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
/// <seealso cref="PDF"/>
public double Density(double x)
{
return DensityImpl(Scale, Shape, Truncation, x);
}
/// <summary>
/// Computes the log probability density of the distribution (lnPDF) at x, i.e. ln(∂P(X ≤ x)/∂x).
/// </summary>
/// <param name="x">The location at which to compute the log density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
/// <seealso cref="PDFLn"/>
public double DensityLn(double x)
{
return DensityLnImpl(Scale, Shape, Truncation, x);
}
/// <summary>
/// Computes the cumulative distribution (CDF) of the distribution at x, i.e. P(X ≤ x).
/// </summary>
/// <param name="x">The location at which to compute the cumulative distribution function.</param>
/// <returns>the cumulative distribution at location <paramref name="x"/>.</returns>
/// <seealso cref="CDF"/>
public double CumulativeDistribution(double x)
{
return CumulativeDistributionImpl(Scale, Shape, Truncation, x);
}
/// <summary>
/// Computes the inverse cumulative distribution (CDF) of the distribution at p, i.e. solving for P(X ≤ x) = p.
/// </summary>
/// <param name="p">The location at which to compute the inverse cumulative distribution function.</param>
/// <returns>the inverse cumulative distribution at location <paramref name="p"/>.</returns>
public double InvCDF(double p)
{
return InvCDFUncheckedImpl(Scale, Shape, Truncation, p);
}
/// <summary>
/// Computes the inverse cumulative distribution (CDF) of the distribution at p, i.e. solving for P(X ≤ x) = p.
/// </summary>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
/// <param name="p">The location at which to compute the inverse cumulative distribution function.</param>
/// <returns>the inverse cumulative distribution at location <paramref name="p"/>.</returns>
/// <seealso cref="CumulativeDistribution"/>
public static double ICDF(double scale, double shape, double truncation, double p)
{
if (!IsValidParameterSet(scale, shape, truncation))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return InvCDFUncheckedImpl(scale, shape, truncation, p);
}
/// <summary>
/// Computes the probability density of the distribution (PDF) at x, i.e. ∂P(X ≤ x)/∂x.
/// </summary>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
/// <seealso cref="Density"/>
public static double PDF(double scale, double shape, double truncation, double x)
{
if (!IsValidParameterSet(scale, shape, truncation))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return DensityImpl(scale, shape, truncation, x);
}
/// <summary>
/// Computes the log probability density of the distribution (lnPDF) at x, i.e. ln(∂P(X ≤ x)/∂x).
/// </summary>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
/// <param name="x">The location at which to compute the log density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
/// <seealso cref="DensityLn"/>
public static double PDFLn(double scale, double shape, double truncation, double x)
{
if (!IsValidParameterSet(scale, shape, truncation))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return DensityLnImpl(scale, shape, truncation, x);
}
/// <summary>
/// Computes the cumulative distribution (CDF) of the distribution at x, i.e. P(X ≤ x).
/// </summary>
/// <param name="scale">The scale (xm) of the distribution. Range: xm > 0.</param>
/// <param name="shape">The shape (α) of the distribution. Range: α > 0.</param>
/// <param name="truncation">The truncation (T) of the distribution. Range: T > xm.</param>
/// <param name="x">The location at which to compute the cumulative distribution function.</param>
/// <returns>the cumulative distribution at location <paramref name="x"/>.</returns>
/// <seealso cref="CumulativeDistribution"/>
public static double CDF(double scale, double shape, double truncation, double x)
{
if (!IsValidParameterSet(scale, shape, truncation))
{
throw new ArgumentException(Resources.InvalidDistributionParameters);
}
return CumulativeDistributionImpl(scale, shape, truncation, x);
}
internal static double DensityImpl(double scale, double shape, double truncation, double x)
{
if (x < scale || x > truncation)
return 0;
else
return (shape * Math.Pow(scale, shape) * Math.Pow(x, -shape - 1)) / (1 - Math.Pow(scale / truncation, shape));
}
internal static double DensityLnImpl(double scale, double shape, double truncation, double x)
{
return Math.Log(DensityImpl(scale, shape, truncation, x));
}
internal static double CumulativeDistributionImpl(double scale, double shape, double truncation, double x)
{
if (x <= scale)
return 0;
else if (x >= truncation)
return 1;
else
return (1 - Math.Pow(scale, shape) * Math.Pow(x, -shape)) / (1 - Math.Pow(scale / truncation, shape));
}
internal static double InvCDFUncheckedImpl(double scale, double shape, double truncation, double p)
{
var numerator = p * Math.Pow(truncation, shape) - p * Math.Pow(scale, shape) - Math.Pow(truncation, shape);
var denominator = Math.Pow(truncation, shape) * Math.Pow(scale, shape);
return Math.Pow(-numerator / denominator, -(1 / shape));
}
}
}

9
src/Numerics/Precision.Comparison.cs

@ -674,5 +674,14 @@ namespace MathNet.Numerics
return CompareToNumbersBetween(a, b, maxNumbersBetween) < 0; return CompareToNumbersBetween(a, b, maxNumbersBetween) < 0;
} }
/// <summary>
/// Checks if a given double values is finite, i.e. neither NaN nor inifnity
/// </summary>
/// <param name="value">The value to be checked fo finitenes.</param>
public static bool IsFinite(this double value)
{
return !double.IsNaN(value) && !double.IsInfinity(value);
}
} }
} }

11
src/Numerics/Properties/Resources.Designer.cs

@ -853,6 +853,17 @@ namespace MathNet.Numerics.Properties
} }
} }
/// <summary>
/// Looks up a localized string similar to: Numerical estimation of the statistic has failed.
/// </summary>
public static string NumericalEstimationFailed
{
get
{
return ResourceManager.GetString("NumericalEstimationFailed", resourceCulture);
}
}
/// <summary> /// <summary>
/// Looks up a localized string similar to The two arguments can&apos;t be compared (maybe they are part of a partial ordering?). /// Looks up a localized string similar to The two arguments can&apos;t be compared (maybe they are part of a partial ordering?).
/// </summary> /// </summary>

3
src/Numerics/Properties/Resources.resx

@ -448,4 +448,7 @@
<data name="DegreesOfFreedomMustBeLessThanSampleSize" xml:space="preserve"> <data name="DegreesOfFreedomMustBeLessThanSampleSize" xml:space="preserve">
<value>The sample size must be larger than the given degrees of freedom.</value> <value>The sample size must be larger than the given degrees of freedom.</value>
</data> </data>
<data name="NumericalEstimationFailed" xml:space="preserve">
<value>Numerical estimation of the statistic has failed. The used solver did not succeed in finding a root.</value>
</data>
</root> </root>
Loading…
Cancel
Save