ufkc7e3d
/
mathnet-numerics
forked from tsai/mathnet-numerics
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity
Browse Source

Distributions: log-normal InvCDF, tweak prototype

optimization-1
Christoph Ruegg 13 years ago
parent
f1dbedd766
commit
81ba0967b0
3 changed files with 87 additions and 11 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Split View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 50
      src/Numerics/Distributions/LogNormal.cs
  2. 14
      src/Numerics/Distributions/Normal.cs
  3. 34
      src/UnitTests/DistributionTests/Continuous/LogNormalTests.cs

50
src/Numerics/Distributions/LogNormal.cs
View File

@ -259,6 +259,7 @@ namespace MathNet.Numerics.Distributions
/// </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)
{
if (x < 0.0)
@ -275,6 +276,7 @@ namespace MathNet.Numerics.Distributions
/// </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)
{
if (x < 0.0)
@ -291,11 +293,24 @@ namespace MathNet.Numerics.Distributions
/// </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 x < 0.0
? 0.0
: 0.5*(1.0 + SpecialFunctions.Erf((Math.Log(x) - _mu)/(_sigma*Constants.Sqrt2)));
return x < 0.0 ? 0.0
: 0.5*SpecialFunctions.Erfc((_mu - Math.Log(x))/(_sigma*Constants.Sqrt2));
}
/// <summary>
/// Computes the inverse of the cumulative distribution function (InvCDF) for the distribution
/// at the given probability. This is also known as the 'quantile function'.
/// </summary>
/// <param name="p">The location at which to compute the inverse cumulative density.</param>
/// <returns>the inverse cumulative density at <paramref name="p"/>.</returns>
/// <seealso cref="InvCDF"/>
public double InverseCumulativeDistribution(double p)
{
return p <= 0.0 ? 0.0 : p >= 1.0 ? double.PositiveInfinity
: Math.Exp(_mu - _sigma*Constants.Sqrt2*SpecialFunctions.ErfcInv(2.0*p));
}
/// <summary>
@ -323,9 +338,10 @@ namespace MathNet.Numerics.Distributions
/// <param name="mu">The log-scale (μ) of the distribution.</param>
/// <param name="sigma">The shape (σ) of the distribution. Range: σ ≥ 0.</param>
/// <returns>the density at <paramref name="x"/>.</returns>
/// <seealso cref="Density"/>
public static double PDF(double mu, double sigma, double x)
{
if (sigma < 0.0) throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
if (sigma < 0.0) throw new ArgumentOutOfRangeException("sigma", Resources.InvalidDistributionParameters);
if (x < 0.0)
{
@ -343,9 +359,10 @@ namespace MathNet.Numerics.Distributions
/// <param name="mu">The log-scale (μ) of the distribution.</param>
/// <param name="sigma">The shape (σ) of the distribution. Range: σ ≥ 0.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
/// <seealso cref="DensityLn"/>
public static double PDFLn(double mu, double sigma, double x)
{
if (sigma < 0.0) throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
if (sigma < 0.0) throw new ArgumentOutOfRangeException("sigma", Resources.InvalidDistributionParameters);
if (x < 0.0)
{
@ -363,15 +380,32 @@ namespace MathNet.Numerics.Distributions
/// <param name="mu">The log-scale (μ) of the distribution.</param>
/// <param name="sigma">The shape (σ) of the distribution. Range: σ ≥ 0.</param>
/// <returns>the cumulative distribution at location <paramref name="x"/>.</returns>
/// <seealso cref="CumulativeDistribution"/>
public static double CDF(double mu, double sigma, double x)
{
if (sigma < 0.0) throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
if (sigma < 0.0) throw new ArgumentOutOfRangeException("sigma", Resources.InvalidDistributionParameters);
return x < 0.0
? 0.0
return x < 0.0 ? 0.0
: 0.5*(1.0 + SpecialFunctions.Erf((Math.Log(x) - mu)/(sigma*Constants.Sqrt2)));
}
/// <summary>
/// Computes the inverse of the cumulative distribution function (InvCDF) for the distribution
/// at the given probability. This is also known as the 'quantile function'.
/// </summary>
/// <param name="p">The location at which to compute the inverse cumulative density.</param>
/// <param name="mu">The log-scale (μ) of the distribution.</param>
/// <param name="sigma">The shape (σ) of the distribution. Range: σ ≥ 0.</param>
/// <returns>the inverse cumulative density at <paramref name="p"/>.</returns>
/// <seealso cref="InverseCumulativeDistribution"/>
public static double InvCDF(double mu, double sigma, double p)
{
if (sigma < 0.0) throw new ArgumentOutOfRangeException("sigma", Resources.InvalidDistributionParameters);
return p <= 0.0 ? 0.0 : p >= 1.0 ? double.PositiveInfinity
: Math.Exp(mu - sigma*Constants.Sqrt2*SpecialFunctions.ErfcInv(2.0*p));
}
/// <summary>
/// Generates a sample from the log-normal distribution using the <i>Box-Muller</i> algorithm.
/// </summary>

14
src/Numerics/Distributions/Normal.cs
View File

@ -279,6 +279,7 @@ namespace MathNet.Numerics.Distributions
/// </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)
{
var d = (x - _mean)/_stdDev;
@ -290,6 +291,7 @@ namespace MathNet.Numerics.Distributions
/// </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)
{
var d = (x - _mean)/_stdDev;
@ -301,17 +303,19 @@ namespace MathNet.Numerics.Distributions
/// </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 0.5*(1.0 + SpecialFunctions.Erf((x - _mean)/(_stdDev*Constants.Sqrt2)));
return 0.5*SpecialFunctions.Erfc((_mean - x)/(_stdDev*Constants.Sqrt2));
}
/// <summary>
/// Computes the inverse of the cumulative distribution function (InvCDF) for the distribution
/// at the given probability.
/// at the given probability. This is also known as the 'quantile function'.
/// </summary>
/// <param name="p">The location at which to compute the inverse cumulative density.</param>
/// <returns>the inverse cumulative density at <paramref name="p"/>.</returns>
/// <seealso cref="InvCDF"/>
public double InverseCumulativeDistribution(double p)
{
return _mean - (_stdDev*Constants.Sqrt2*SpecialFunctions.ErfcInv(2.0*p));
@ -368,6 +372,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="stddev">The standard deviation (σ) of the normal 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 mean, double stddev, double x)
{
if (stddev < 0.0) throw new ArgumentOutOfRangeException("stddev", Resources.InvalidDistributionParameters);
@ -383,6 +388,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="stddev">The standard deviation (σ) of the normal distribution. Range: σ ≥ 0.</param>
/// <param name="x">The location at which to compute the density.</param>
/// <returns>the log density at <paramref name="x"/>.</returns>
/// <seealso cref="DensityLn"/>
public static double PDFLn(double mean, double stddev, double x)
{
if (stddev < 0.0) throw new ArgumentOutOfRangeException("stddev", Resources.InvalidDistributionParameters);
@ -398,6 +404,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="mean">The mean (μ) of the normal distribution.</param>
/// <param name="stddev">The standard deviation (σ) of the normal distribution. Range: σ ≥ 0.</param>
/// <returns>the cumulative distribution at location <paramref name="x"/>.</returns>
/// <seealso cref="CumulativeDistribution"/>
public static double CDF(double mean, double stddev, double x)
{
if (stddev < 0.0) throw new ArgumentOutOfRangeException("stddev", Resources.InvalidDistributionParameters);
@ -407,12 +414,13 @@ namespace MathNet.Numerics.Distributions
/// <summary>
/// Computes the inverse of the cumulative distribution function (InvCDF) for the distribution
/// at the given probability.
/// at the given probability. This is also known as the 'quantile function'.
/// </summary>
/// <param name="p">The location at which to compute the inverse cumulative density.</param>
/// <param name="mean">The mean (μ) of the normal distribution.</param>
/// <param name="stddev">The standard deviation (σ) of the normal distribution. Range: σ ≥ 0.</param>
/// <returns>the inverse cumulative density at <paramref name="p"/>.</returns>
/// <seealso cref="InverseCumulativeDistribution"/>
public static double InvCDF(double mean, double stddev, double p)
{
if (stddev < 0.0) throw new ArgumentOutOfRangeException("stddev", Resources.InvalidDistributionParameters);

34
src/UnitTests/DistributionTests/Continuous/LogNormalTests.cs
View File

@ -523,6 +523,40 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
AssertHelpers.AlmostEqual(f, LogNormal.CDF(mu, sigma, x), 8);
}
/// <summary>
/// Validate inverse cumulative distribution.
/// </summary>
/// <param name="mu">Mu parameter.</param>
/// <param name="sigma">Sigma value.</param>
/// <param name="x">Input X value.</param>
/// <param name="f">Expected value.</param>
[TestCase(-0.100000, 0.100000, 0.500000, 0.0000000015011556178148777579869633555518882664666520593658)]
[TestCase(-0.100000, 0.100000, 0.800000, 0.10908001076375810900224507908874442583171381706127)]
[TestCase(-0.100000, 1.500000, 0.100000, 0.070999149762464508991968731574953594549291668468349)]
[TestCase(-0.100000, 1.500000, 0.500000, 0.34626224992888089297789445771047690175505847991946)]
[TestCase(-0.100000, 1.500000, 0.800000, 0.46728530589487698517090261668589508746353129242404)]
[TestCase(-0.100000, 2.500000, 0.100000, 0.18914969879695093477606645992572208111152994999076)]
[TestCase(-0.100000, 2.500000, 0.500000, 0.40622798321378106125020505907901206714868922279347)]
[TestCase(-0.100000, 2.500000, 0.800000, 0.48035707589956665425068652807400957345208517749893)]
[TestCase(1.500000, 1.500000, 0.100000, 0.005621455876973168709588070988239748831823850202953)]
[TestCase(1.500000, 1.500000, 0.500000, 0.07185716187918271235246980951571040808235628115265)]
[TestCase(1.500000, 1.500000, 0.800000, 0.12532699044614938400496547188720940854423187977236)]
[TestCase(1.500000, 2.500000, 0.100000, 0.064125647996943514411570834861724406903677144126117)]
[TestCase(1.500000, 2.500000, 0.500000, 0.19017302281590810871719754032332631806011441356498)]
[TestCase(1.500000, 2.500000, 0.800000, 0.24533064397555500690927047163085419096928289095201)]
[TestCase(2.500000, 1.500000, 0.100000, 0.00068304052220788502001572635016579586444611070077399)]
[TestCase(2.500000, 1.500000, 0.500000, 0.016636862816580533038130583128179878924863968664206)]
[TestCase(2.500000, 1.500000, 0.800000, 0.034729001282904174941366974418836262996834852343018)]
[TestCase(2.500000, 2.500000, 0.100000, 0.027363708266690978870139978537188410215717307180775)]
[TestCase(2.500000, 2.500000, 0.500000, 0.10075543423327634536450625420610429181921642201567)]
[TestCase(2.500000, 2.500000, 0.800000, 0.13802019192453118732001307556787218421918336849121)]
public void ValidateInverseCumulativeDistribution(double mu, double sigma, double x, double f)
{
var n = new LogNormal(mu, sigma);
AssertHelpers.AlmostEqual(x, n.InverseCumulativeDistribution(f), 8);
AssertHelpers.AlmostEqual(x, LogNormal.InvCDF(mu, sigma, f), 8);
}
/// <summary>
/// Can estimate distribution parameters.
/// </summary>

Write Preview
Loading…
Cancel
Save