Browse Source

Distributions: consistent sample methods for continuous distributions #32

la-knuth
Christoph Ruegg 14 years ago
parent
commit
6c32ba172b
  1. 140
      src/Numerics/Distributions/Continuous/Beta.cs
  2. 136
      src/Numerics/Distributions/Continuous/Cauchy.cs
  3. 113
      src/Numerics/Distributions/Continuous/Chi.cs
  4. 91
      src/Numerics/Distributions/Continuous/ChiSquare.cs
  5. 79
      src/Numerics/Distributions/Continuous/ContinuousUniform.cs
  6. 174
      src/Numerics/Distributions/Continuous/Erlang.cs
  7. 115
      src/Numerics/Distributions/Continuous/Exponential.cs
  8. 111
      src/Numerics/Distributions/Continuous/FisherSnedecor.cs
  9. 213
      src/Numerics/Distributions/Continuous/Gamma.cs
  10. 126
      src/Numerics/Distributions/Continuous/InverseGamma.cs
  11. 135
      src/Numerics/Distributions/Continuous/Laplace.cs
  12. 77
      src/Numerics/Distributions/Continuous/LogNormal.cs
  13. 179
      src/Numerics/Distributions/Continuous/Normal.cs
  14. 123
      src/Numerics/Distributions/Continuous/Pareto.cs
  15. 118
      src/Numerics/Distributions/Continuous/Rayleigh.cs
  16. 164
      src/Numerics/Distributions/Continuous/Stable.cs
  17. 131
      src/Numerics/Distributions/Continuous/StudentT.cs
  18. 103
      src/Numerics/Distributions/Continuous/Weibull.cs
  19. 2
      src/Numerics/Distributions/Multivariate/MatrixNormal.cs

140
src/Numerics/Distributions/Continuous/Beta.cs

@ -51,17 +51,17 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Beta shape parameter a. /// Beta shape parameter a.
/// </summary> /// </summary>
private double _shapeA; double _shapeA;
/// <summary> /// <summary>
/// Beta shape parameter b. /// Beta shape parameter b.
/// </summary> /// </summary>
private double _shapeB; double _shapeB;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the Beta class. /// Initializes a new instance of the Beta class.
@ -90,7 +90,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="a">The a shape parameter of the Beta distribution.</param> /// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param> /// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double a, double b) static bool IsValidParameterSet(double a, double b)
{ {
if (a < 0.0 || b < 0.0 || Double.IsNaN(a) || Double.IsNaN(b)) if (a < 0.0 || b < 0.0 || Double.IsNaN(a) || Double.IsNaN(b))
{ {
@ -106,7 +106,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="a">The a shape parameter of the Beta distribution.</param> /// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param> /// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double a, double b) void SetParameters(double a, double b)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b)) if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b))
{ {
@ -122,15 +122,8 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double A public double A
{ {
get get { return _shapeA; }
{ set { SetParameters(value, _shapeB); }
return _shapeA;
}
set
{
SetParameters(value, _shapeB);
}
} }
/// <summary> /// <summary>
@ -138,15 +131,8 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double B public double B
{ {
get get { return _shapeB; }
{ set { SetParameters(_shapeA, value); }
return _shapeB;
}
set
{
SetParameters(_shapeA, value);
}
} }
#region IDistribution implementation #region IDistribution implementation
@ -156,10 +142,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -218,10 +201,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get get { return (_shapeA * _shapeB) / ((_shapeA + _shapeB) * (_shapeA + _shapeB) * (_shapeA + _shapeB + 1.0)); }
{
return (_shapeA * _shapeB) / ((_shapeA + _shapeB) * (_shapeA + _shapeB) * (_shapeA + _shapeB + 1.0));
}
} }
/// <summary> /// <summary>
@ -229,10 +209,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return Math.Sqrt((_shapeA * _shapeB) / ((_shapeA + _shapeB) * (_shapeA + _shapeB) * (_shapeA + _shapeB + 1.0))); }
{
return Math.Sqrt((_shapeA * _shapeB) / ((_shapeA + _shapeB) * (_shapeA + _shapeB) * (_shapeA + _shapeB + 1.0)));
}
} }
/// <summary> /// <summary>
@ -258,9 +235,9 @@ namespace MathNet.Numerics.Distributions
} }
return SpecialFunctions.BetaLn(_shapeA, _shapeB) return SpecialFunctions.BetaLn(_shapeA, _shapeB)
- ((_shapeA - 1.0) * SpecialFunctions.DiGamma(_shapeA)) - ((_shapeA - 1.0) * SpecialFunctions.DiGamma(_shapeA))
- ((_shapeB - 1.0) * SpecialFunctions.DiGamma(_shapeB)) - ((_shapeB - 1.0) * SpecialFunctions.DiGamma(_shapeB))
+ ((_shapeA + _shapeB - 2.0) * SpecialFunctions.DiGamma(_shapeA + _shapeB)); + ((_shapeA + _shapeB - 2.0) * SpecialFunctions.DiGamma(_shapeA + _shapeB));
} }
} }
@ -302,7 +279,7 @@ namespace MathNet.Numerics.Distributions
} }
return 2.0 * (_shapeB - _shapeA) * Math.Sqrt(_shapeA + _shapeB + 1.0) return 2.0 * (_shapeB - _shapeA) * Math.Sqrt(_shapeA + _shapeB + 1.0)
/ ((_shapeA + _shapeB + 2.0) * Math.Sqrt(_shapeA * _shapeB)); / ((_shapeA + _shapeB + 2.0) * Math.Sqrt(_shapeA * _shapeB));
} }
} }
@ -361,10 +338,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -372,10 +346,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -383,10 +354,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return 1.0; }
{
return 1.0;
}
} }
/// <summary> /// <summary>
@ -515,27 +483,27 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (x >= 1.0) if (x >= 1.0)
{ {
return 1.0; return 1.0;
} }
if (Double.IsPositiveInfinity(_shapeA) && Double.IsPositiveInfinity(_shapeB)) if (Double.IsPositiveInfinity(_shapeA) && Double.IsPositiveInfinity(_shapeB))
{ {
return x < 0.5 ? 0.0 : 1.0; return x < 0.5 ? 0.0 : 1.0;
} }
if (Double.IsPositiveInfinity(_shapeA)) if (Double.IsPositiveInfinity(_shapeA))
{ {
return x < 1.0 ? 0.0 : 1.0; return x < 1.0 ? 0.0 : 1.0;
} }
if (Double.IsPositiveInfinity(_shapeB)) if (Double.IsPositiveInfinity(_shapeB))
{ {
return x >= 0.0 ? 1.0 : 0.0; return x >= 0.0 ? 1.0 : 0.0;
} }
if (_shapeA == 0.0 && _shapeB == 0.0) if (_shapeA == 0.0 && _shapeB == 0.0)
{ {
if (x >= 0.0 && x < 1.0) if (x >= 0.0 && x < 1.0)
@ -545,32 +513,48 @@ namespace MathNet.Numerics.Distributions
return 1.0; return 1.0;
} }
if (_shapeA == 0.0) if (_shapeA == 0.0)
{ {
return 1.0; return 1.0;
} }
if (_shapeB == 0.0) if (_shapeB == 0.0)
{ {
return x >= 1.0 ? 1.0 : 0.0; return x >= 1.0 ? 1.0 : 0.0;
} }
if (_shapeA == 1.0 && _shapeB == 1.0) if (_shapeA == 1.0 && _shapeB == 1.0)
{ {
return x; return x;
} }
return SpecialFunctions.BetaRegularized(_shapeA, _shapeB, x); return SpecialFunctions.BetaRegularized(_shapeA, _shapeB, x);
} }
#endregion
/// <summary>
/// Samples Beta distributed random variables by sampling two Gamma variables and normalizing.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="a">The A shape parameter.</param>
/// <param name="b">The B shape parameter.</param>
/// <returns>a random number from the Beta distribution.</returns>
internal static double SampleUnchecked(Random rnd, double a, double b)
{
var x = Gamma.SampleUnchecked(rnd, a, 1.0);
var y = Gamma.SampleUnchecked(rnd, b, 1.0);
return x / (x + y);
}
/// <summary> /// <summary>
/// Generates a sample from the Beta distribution. /// Generates a sample from the Beta distribution.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return SampleBeta(RandomSource, _shapeA, _shapeB); return SampleUnchecked(RandomSource, _shapeA, _shapeB);
} }
/// <summary> /// <summary>
@ -581,37 +565,35 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return SampleBeta(RandomSource, _shapeA, _shapeB); yield return SampleUnchecked(RandomSource, _shapeA, _shapeB);
} }
} }
#endregion
/// <summary> /// <summary>
/// Generates a sample from the normal distribution using the <i>Box-Muller</i> algorithm. /// Generates a sample from the distribution.
/// </summary> /// </summary>
/// <param name="rng">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="a">The a shape parameter of the Beta distribution.</param> /// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param> /// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public static double Sample(Random rng, double a, double b) public static double Sample(Random rnd, double a, double b)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b)) if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b))
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return SampleBeta(rng, a, b); return SampleUnchecked(rnd, a, b);
} }
/// <summary> /// <summary>
/// Generates a sequence of samples from the normal distribution using the <i>Box-Muller</i> algorithm. /// Generates a sequence of samples from the distribution.
/// </summary> /// </summary>
/// <param name="rng">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="a">The a shape parameter of the Beta distribution.</param> /// <param name="a">The a shape parameter of the Beta distribution.</param>
/// <param name="b">The b shape parameter of the Beta distribution.</param> /// <param name="b">The b shape parameter of the Beta distribution.</param>
/// <returns>a sequence of samples from the distribution.</returns> /// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rng, double a, double b) public static IEnumerable<double> Samples(Random rnd, double a, double b)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b)) if (Control.CheckDistributionParameters && !IsValidParameterSet(a, b))
{ {
@ -620,22 +602,8 @@ namespace MathNet.Numerics.Distributions
while (true) while (true)
{ {
yield return SampleBeta(rng, a, b); yield return SampleUnchecked(rnd, a, b);
} }
} }
/// <summary>
/// Samples Beta distributed random variables by sampling two Gamma variables and normalizing.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="a">The A shape parameter.</param>
/// <param name="b">The B shape parameter.</param>
/// <returns>a random number from the Beta distribution.</returns>
internal static double SampleBeta(Random rnd, double a, double b)
{
var x = Gamma.SampleGamma(rnd, a, 1.0);
var y = Gamma.SampleGamma(rnd, b, 1.0);
return x / (x + y);
}
} }
} }

136
src/Numerics/Distributions/Continuous/Cauchy.cs

@ -44,17 +44,18 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The scale of the Cauchy distribution. /// The scale of the Cauchy distribution.
/// </summary> /// </summary>
private double _scale; double _scale;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="Cauchy"/> class with the location parameter set to 0 and the scale parameter set to 1 /// Initializes a new instance of the <see cref="Cauchy"/> class with the location parameter set to 0 and the scale parameter set to 1
/// </summary> /// </summary>
public Cauchy() : this(0, 1) public Cauchy()
: this(0, 1)
{ {
} }
@ -82,7 +83,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="location">Location parameter.</param> /// <param name="location">Location parameter.</param>
/// <param name="scale">Scale parameter. Must be greater than 0.</param> /// <param name="scale">Scale parameter. Must be greater than 0.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double location, double scale) void SetParameters(double location, double scale)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale)) if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale))
{ {
@ -99,7 +100,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="location">Location parameter.</param> /// <param name="location">Location parameter.</param>
/// <param name="scale">Scale parameter. Must be greater than 0.</param> /// <param name="scale">Scale parameter. Must be greater than 0.</param>
/// <returns>True when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns>True when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double location, double scale) static bool IsValidParameterSet(double location, double scale)
{ {
if (scale <= 0) if (scale <= 0)
{ {
@ -119,15 +120,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Location public double Location
{ {
get get { return Median; }
{
return Median;
}
set set { SetParameters(value, _scale); }
{
SetParameters(value, _scale);
}
} }
/// <summary> /// <summary>
@ -135,15 +130,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return _scale; }
{
return _scale;
}
set set { SetParameters(Median, value); }
{
SetParameters(Median, value);
}
} }
/// <summary> /// <summary>
@ -162,11 +151,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
if (value == null) if (value == null)
@ -183,10 +168,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -194,10 +176,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -205,10 +184,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -216,10 +192,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return Math.Log(4.0 * Constants.Pi * _scale); }
{
return Math.Log(4.0 * Constants.Pi * _scale);
}
} }
/// <summary> /// <summary>
@ -227,10 +200,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Skewness public double Skewness
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -252,30 +222,20 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return Median; }
{
return Median;
}
} }
/// <summary> /// <summary>
/// Gets the median of the distribution. /// Gets the median of the distribution.
/// </summary> /// </summary>
public double Median public double Median { get; private set; }
{
get;
private set;
}
/// <summary> /// <summary>
/// Gets the minimum of the distribution. /// Gets the minimum of the distribution.
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return Double.NegativeInfinity; }
{
return Double.NegativeInfinity;
}
} }
/// <summary> /// <summary>
@ -283,10 +243,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -309,13 +266,28 @@ namespace MathNet.Numerics.Distributions
return -Math.Log(Constants.Pi * _scale * (1.0 + (((x - Median) / _scale) * ((x - Median) / _scale)))); return -Math.Log(Constants.Pi * _scale * (1.0 + (((x - Median) / _scale) * ((x - Median) / _scale))));
} }
#endregion
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="location">The location shape parameter.</param>
/// <param name="scale">The scale parameter.</param>
/// <returns>a random number from the distribution.</returns>
internal static double SampleUnchecked(Random rnd, double location, double scale)
{
var u = rnd.NextDouble();
return location + (scale * Math.Tan(Constants.Pi * (u - 0.5)));
}
/// <summary> /// <summary>
/// Draws a random sample from the distribution. /// Draws a random sample from the distribution.
/// </summary> /// </summary>
/// <returns>A random number from this distribution.</returns> /// <returns>A random number from this distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, Median, _scale); return SampleUnchecked(RandomSource, Median, _scale);
} }
/// <summary> /// <summary>
@ -326,23 +298,45 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource, Median, _scale); yield return SampleUnchecked(RandomSource, Median, _scale);
} }
} }
#endregion /// <summary>
/// Generates a sample from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="location">The location shape parameter.</param>
/// <param name="scale">The scale parameter.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rnd, double location, double scale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, location, scale);
}
/// <summary> /// <summary>
/// Samples the distribution. /// Generates a sequence of samples from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="location">The location shape parameter.</param> /// <param name="location">The location shape parameter.</param>
/// <param name="scale">The scale parameter.</param> /// <param name="scale">The scale parameter.</param>
/// <returns>a random number from the distribution.</returns> /// <returns>a sequence of samples from the distribution.</returns>
private static double DoSample(Random rnd, double location, double scale) public static IEnumerable<double> Samples(Random rnd, double location, double scale)
{ {
var u = rnd.NextDouble(); if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale))
return location + (scale * Math.Tan(Constants.Pi * (u - 0.5))); {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleUnchecked(rnd, location, scale);
}
} }
} }
} }

113
src/Numerics/Distributions/Continuous/Chi.cs

@ -47,12 +47,12 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Keeps track of the degrees of freedom for the Chi distribution. /// Keeps track of the degrees of freedom for the Chi distribution.
/// </summary> /// </summary>
private double _dof; double _dof;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="Chi"/> class. /// Initializes a new instance of the <see cref="Chi"/> class.
@ -71,7 +71,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="dof">The degrees of freedom for the Chi distribution.</param> /// <param name="dof">The degrees of freedom for the Chi distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double dof) void SetParameters(double dof)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(dof)) if (Control.CheckDistributionParameters && !IsValidParameterSet(dof))
{ {
@ -86,7 +86,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="dof">The degrees of freedom for the Chi distribution.</param> /// <param name="dof">The degrees of freedom for the Chi distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double dof) static bool IsValidParameterSet(double dof)
{ {
if (dof <= 0 || Double.IsNaN(dof)) if (dof <= 0 || Double.IsNaN(dof))
{ {
@ -101,15 +101,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double DegreesOfFreedom public double DegreesOfFreedom
{ {
get get { return _dof; }
{
return _dof;
}
set set { SetParameters(value); }
{
SetParameters(value);
}
} }
/// <summary> /// <summary>
@ -128,10 +122,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -149,10 +140,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get get { return Math.Sqrt(2) * (SpecialFunctions.Gamma((_dof + 1.0) / 2.0) / SpecialFunctions.Gamma(_dof / 2.0)); }
{
return Math.Sqrt(2) * (SpecialFunctions.Gamma((_dof + 1.0) / 2.0) / SpecialFunctions.Gamma(_dof / 2.0));
}
} }
/// <summary> /// <summary>
@ -160,10 +148,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get get { return _dof - (Mean * Mean); }
{
return _dof - (Mean * Mean);
}
} }
/// <summary> /// <summary>
@ -171,10 +156,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return Math.Sqrt(Variance); }
{
return Math.Sqrt(Variance);
}
} }
/// <summary> /// <summary>
@ -182,10 +164,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return SpecialFunctions.GammaLn(_dof / 2.0) + ((_dof - Math.Log(2) - ((_dof - 1.0) * SpecialFunctions.DiGamma(_dof / 2.0))) / 2.0); }
{
return SpecialFunctions.GammaLn(_dof / 2.0) + ((_dof - Math.Log(2) - ((_dof - 1.0) * SpecialFunctions.DiGamma(_dof / 2.0))) / 2.0);
}
} }
/// <summary> /// <summary>
@ -235,10 +214,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -246,10 +222,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -257,10 +230,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -283,13 +253,31 @@ namespace MathNet.Numerics.Distributions
return ((1.0 - (_dof / 2.0)) * Math.Log(2.0)) + ((_dof - 1.0) * Math.Log(x)) - (x * x / 2.0) - SpecialFunctions.GammaLn(_dof / 2.0); return ((1.0 - (_dof / 2.0)) * Math.Log(2.0)) + ((_dof - 1.0) * Math.Log(x)) - (x * x / 2.0) - SpecialFunctions.GammaLn(_dof / 2.0);
} }
#endregion
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <returns>a random number from the distribution.</returns>
internal static double SampleUnchecked(Random rnd, int dof)
{
double sum = 0;
for (var i = 0; i < dof; i++)
{
sum += Math.Pow(Normal.Sample(rnd, 0.0, 1.0), 2);
}
return Math.Sqrt(sum);
}
/// <summary> /// <summary>
/// Generates a sample from the Chi distribution. /// Generates a sample from the Chi distribution.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource); return SampleUnchecked(RandomSource, (int)_dof);
} }
/// <summary> /// <summary>
@ -298,29 +286,44 @@ namespace MathNet.Numerics.Distributions
/// <returns>a sequence of samples from the distribution.</returns> /// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples() public IEnumerable<double> Samples()
{ {
var dof = (int)_dof;
while (true) while (true)
{ {
yield return DoSample(RandomSource); yield return SampleUnchecked(RandomSource, dof);
} }
} }
#endregion /// <summary>
/// Generates a sample from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rnd, int dof)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(dof))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, dof);
}
/// <summary> /// <summary>
/// Samples the distribution. /// Generates a sequence of samples from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <returns>a random number from the distribution.</returns> /// <returns>a sequence of samples from the distribution.</returns>
private double DoSample(Random rnd) public static IEnumerable<double> Samples(Random rnd, int dof)
{ {
double sum = 0; if (Control.CheckDistributionParameters && !IsValidParameterSet(dof))
var n = (int)_dof;
for (var i = 0; i < n; i++)
{ {
sum += Math.Pow(Normal.Sample(rnd, 0.0, 1.0), 2); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return Math.Sqrt(sum); while (true)
{
yield return SampleUnchecked(rnd, dof);
}
} }
} }
} }

91
src/Numerics/Distributions/Continuous/ChiSquare.cs

@ -49,7 +49,7 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="ChiSquare"/> class. /// Initializes a new instance of the <see cref="ChiSquare"/> class.
@ -68,7 +68,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="dof">The degrees of freedom for the <c>ChiSquare</c> distribution.</param> /// <param name="dof">The degrees of freedom for the <c>ChiSquare</c> distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double dof) void SetParameters(double dof)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(dof)) if (Control.CheckDistributionParameters && !IsValidParameterSet(dof))
{ {
@ -83,7 +83,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="dof">The degrees of freedom for the <c>ChiSquare</c> distribution.</param> /// <param name="dof">The degrees of freedom for the <c>ChiSquare</c> distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double dof) static bool IsValidParameterSet(double dof)
{ {
return dof > 0 && !Double.IsNaN(dof); return dof > 0 && !Double.IsNaN(dof);
} }
@ -93,15 +93,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double DegreesOfFreedom public double DegreesOfFreedom
{ {
get get { return Mean; }
{
return Mean;
}
set set { SetParameters(value); }
{
SetParameters(value);
}
} }
/// <summary> /// <summary>
@ -120,11 +114,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
if (value == null) if (value == null)
@ -139,11 +129,7 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Gets the mean of the distribution. /// Gets the mean of the distribution.
/// </summary> /// </summary>
public double Mean public double Mean { get; private set; }
{
get;
private set;
}
/// <summary> /// <summary>
/// Gets the variance of the distribution. /// Gets the variance of the distribution.
@ -243,13 +229,39 @@ namespace MathNet.Numerics.Distributions
return (-x / 2.0) + (((Mean / 2.0) - 1.0) * Math.Log(x)) - ((Mean / 2.0) * Math.Log(2)) - SpecialFunctions.GammaLn(Mean / 2.0); return (-x / 2.0) + (((Mean / 2.0) - 1.0) * Math.Log(x)) - ((Mean / 2.0) * Math.Log(2)) - SpecialFunctions.GammaLn(Mean / 2.0);
} }
#endregion
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="dof">The degrees of freedom.</param>
/// <returns>a random number from the distribution.</returns>
internal static double SampleUnchecked(Random rnd, double dof)
{
//Use the simple method if the dof is an integer anyway
if (Math.Floor(dof) == dof && dof < Int32.MaxValue)
{
double sum = 0;
var n = (int)dof;
for (var i = 0; i < n; i++)
{
sum += Math.Pow(Normal.Sample(rnd, 0.0, 1.0), 2);
}
return sum;
}
//Call the gamma function (see http://en.wikipedia.org/wiki/Gamma_distribution#Specializations
//for a justification)
return Gamma.SampleUnchecked(rnd, dof / 2.0, .5);
}
/// <summary> /// <summary>
/// Generates a sample from the <c>ChiSquare</c> distribution. /// Generates a sample from the <c>ChiSquare</c> distribution.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, Mean); return SampleUnchecked(RandomSource, Mean);
} }
/// <summary> /// <summary>
@ -260,50 +272,43 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource, Mean); yield return SampleUnchecked(RandomSource, Mean);
} }
} }
#endregion
/// <summary> /// <summary>
/// Samples the distribution. /// Generates a sample from the <c>ChiSquare</c> distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="dof">The degrees of freedom.</param> /// <param name="dof">The degrees of freedom.</param>
/// <returns>a random number from the distribution.</returns> /// <returns>a sample from the distribution. </returns>
private static double DoSample(Random rnd, double dof) public static double Sample(Random rnd, double dof)
{ {
//Use the simple method if the dof is an integer anyway if (Control.CheckDistributionParameters && !IsValidParameterSet(dof))
if (Math.Floor(dof) == dof && dof < Int32.MaxValue)
{ {
double sum = 0; throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
var n = (int)dof;
for (var i = 0; i < n; i++)
{
sum += Math.Pow(Normal.Sample(rnd, 0.0, 1.0), 2);
}
return sum;
} }
//Call the gamma function (see http://en.wikipedia.org/wiki/Gamma_distribution#Specializations
//for a justification) return SampleUnchecked(rnd, dof);
return Gamma.Sample(rnd, dof / 2.0, .5);
} }
/// <summary> /// <summary>
/// Generates a sample from the <c>ChiSquare</c> distribution. /// Generates a sequence of samples from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="dof">The degrees of freedom.</param> /// <param name="dof">The degrees of freedom.</param>
/// <returns>a sample from the distribution. </returns> /// <returns>a sample from the distribution. </returns>
public static double Sample(Random rnd, double dof) public static IEnumerable<double> Samples(Random rnd, double dof)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(dof)) if (Control.CheckDistributionParameters && !IsValidParameterSet(dof))
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return DoSample(rnd, dof); while (true)
{
yield return SampleUnchecked(rnd, dof);
}
} }
} }
} }

79
src/Numerics/Distributions/Continuous/ContinuousUniform.cs

@ -48,22 +48,23 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The distribution's lower bound. /// The distribution's lower bound.
/// </summary> /// </summary>
private double _lower; double _lower;
/// <summary> /// <summary>
/// The distribution's upper bound. /// The distribution's upper bound.
/// </summary> /// </summary>
private double _upper; double _upper;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the ContinuousUniform class with lower bound 0 and upper bound 1. /// Initializes a new instance of the ContinuousUniform class with lower bound 0 and upper bound 1.
/// </summary> /// </summary>
public ContinuousUniform() : this(0.0, 1.0) public ContinuousUniform()
: this(0.0, 1.0)
{ {
} }
@ -94,13 +95,13 @@ namespace MathNet.Numerics.Distributions
/// <param name="lower">Lower bound.</param> /// <param name="lower">Lower bound.</param>
/// <param name="upper">Upper bound; must be at least as large as <paramref name="lower"/>.</param> /// <param name="upper">Upper bound; must be at least as large as <paramref name="lower"/>.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double lower, double upper) static bool IsValidParameterSet(double lower, double upper)
{ {
if (upper < lower) if (upper < lower)
{ {
return false; return false;
} }
if (Double.IsNaN(upper) || Double.IsNaN(lower)) if (Double.IsNaN(upper) || Double.IsNaN(lower))
{ {
return false; return false;
@ -115,7 +116,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="lower">Lower bound.</param> /// <param name="lower">Lower bound.</param>
/// <param name="upper">Upper bound; must be at least as large as <paramref name="lower"/>.</param> /// <param name="upper">Upper bound; must be at least as large as <paramref name="lower"/>.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double lower, double upper) void SetParameters(double lower, double upper)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper)) if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper))
{ {
@ -131,15 +132,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Lower public double Lower
{ {
get get { return _lower; }
{
return _lower;
}
set set { SetParameters(value, _upper); }
{
SetParameters(value, _upper);
}
} }
/// <summary> /// <summary>
@ -147,15 +142,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Upper public double Upper
{ {
get get { return _upper; }
{
return _upper;
}
set set { SetParameters(_lower, value); }
{
SetParameters(_lower, value);
}
} }
#region IDistribution Members #region IDistribution Members
@ -165,10 +154,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -221,6 +207,7 @@ namespace MathNet.Numerics.Distributions
{ {
get { return 0.0; } get { return 0.0; }
} }
#endregion #endregion
#region IContinuousDistribution Members #region IContinuousDistribution Members
@ -300,7 +287,7 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (x >= _upper) if (x >= _upper)
{ {
return 1.0; return 1.0;
@ -309,13 +296,27 @@ namespace MathNet.Numerics.Distributions
return (x - _lower) / (_upper - _lower); return (x - _lower) / (_upper - _lower);
} }
#endregion
/// <summary>
/// Generates one sample from the <c>ContinuousUniform</c> distribution without parameter checking.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="lower">The lower bound of the uniform random variable.</param>
/// <param name="upper">The upper bound of the uniform random variable.</param>
/// <returns>a uniformly distributed random number.</returns>
internal static double SampleUnchecked(Random rnd, double lower, double upper)
{
return lower + (rnd.NextDouble() * (upper - lower));
}
/// <summary> /// <summary>
/// Generates a sample from the <c>ContinuousUniform</c> distribution. /// Generates a sample from the <c>ContinuousUniform</c> distribution.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, _lower, _upper); return SampleUnchecked(RandomSource, _lower, _upper);
} }
/// <summary> /// <summary>
@ -326,12 +327,10 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource, _lower, _upper); yield return SampleUnchecked(RandomSource, _lower, _upper);
} }
} }
#endregion
/// <summary> /// <summary>
/// Generates a sample from the <c>ContinuousUniform</c> distribution. /// Generates a sample from the <c>ContinuousUniform</c> distribution.
/// </summary> /// </summary>
@ -346,7 +345,7 @@ namespace MathNet.Numerics.Distributions
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return DoSample(rnd, lower, upper); return SampleUnchecked(rnd, lower, upper);
} }
/// <summary> /// <summary>
@ -365,20 +364,8 @@ namespace MathNet.Numerics.Distributions
while (true) while (true)
{ {
yield return DoSample(rnd, lower, upper); yield return SampleUnchecked(rnd, lower, upper);
} }
} }
/// <summary>
/// Generates one sample from the <c>ContinuousUniform</c> distribution without parameter checking.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="lower">The lower bound of the uniform random variable.</param>
/// <param name="upper">The upper bound of the uniform random variable.</param>
/// <returns>a uniformly distributed random number.</returns>
private static double DoSample(Random rnd, double lower, double upper)
{
return lower + (rnd.NextDouble() * (upper - lower));
}
} }
} }

174
src/Numerics/Distributions/Continuous/Erlang.cs

@ -46,17 +46,17 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Erlang shape parameter. /// Erlang shape parameter.
/// </summary> /// </summary>
private double _shape; double _shape;
/// <summary> /// <summary>
/// Erlang inverse scale parameter. /// Erlang inverse scale parameter.
/// </summary> /// </summary>
private double _invScale; double _invScale;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="Erlang"/> class. /// Initializes a new instance of the <see cref="Erlang"/> class.
@ -102,7 +102,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="shape">The shape of the Erlang distribution.</param> /// <param name="shape">The shape of the Erlang distribution.</param>
/// <param name="invScale">The inverse scale of the Erlang distribution.</param> /// <param name="invScale">The inverse scale of the Erlang distribution.</param>
private void SetParameters(double shape, double invScale) void SetParameters(double shape, double invScale)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale)) if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{ {
@ -119,7 +119,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="shape">The shape of the Erlang distribution.</param> /// <param name="shape">The shape of the Erlang distribution.</param>
/// <param name="invScale">The inverse scale of the Erlang distribution.</param> /// <param name="invScale">The inverse scale of the Erlang distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double shape, double invScale) static bool IsValidParameterSet(double shape, double invScale)
{ {
if (shape < 0.0 || invScale < 0.0 || Double.IsNaN(shape) || Double.IsNaN(invScale)) if (shape < 0.0 || invScale < 0.0 || Double.IsNaN(shape) || Double.IsNaN(invScale))
{ {
@ -134,15 +134,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public int Shape public int Shape
{ {
get get { return (int)_shape; }
{
return (int)_shape;
}
set set { SetParameters(value, _invScale); }
{
SetParameters(value, _invScale);
}
} }
/// <summary> /// <summary>
@ -150,10 +144,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return 1.0 / _invScale; }
{
return 1.0 / _invScale;
}
set set
{ {
@ -173,15 +164,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double InvScale public double InvScale
{ {
get get { return _invScale; }
{
return _invScale;
}
set set { SetParameters(_shape, value); }
{
SetParameters(_shape, value);
}
} }
/// <summary> /// <summary>
@ -200,10 +185,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -227,12 +209,12 @@ namespace MathNet.Numerics.Distributions
{ {
return _shape; return _shape;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return _shape / _invScale; return _shape / _invScale;
} }
} }
@ -248,12 +230,12 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return _shape / (_invScale * _invScale); return _shape / (_invScale * _invScale);
} }
} }
@ -269,12 +251,12 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return Math.Sqrt(_shape) / _invScale; return Math.Sqrt(_shape) / _invScale;
} }
} }
@ -295,7 +277,7 @@ namespace MathNet.Numerics.Distributions
{ {
return Double.NaN; return Double.NaN;
} }
return _shape - Math.Log(_invScale) + SpecialFunctions.GammaLn(_shape) + ((1.0 - _shape) * SpecialFunctions.DiGamma(_shape)); return _shape - Math.Log(_invScale) + SpecialFunctions.GammaLn(_shape) + ((1.0 - _shape) * SpecialFunctions.DiGamma(_shape));
} }
} }
@ -311,12 +293,12 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return 2.0 / Math.Sqrt(_shape); return 2.0 / Math.Sqrt(_shape);
} }
} }
@ -332,12 +314,12 @@ namespace MathNet.Numerics.Distributions
{ {
return x >= _shape ? 1.0 : 0.0; return x >= _shape ? 1.0 : 0.0;
} }
if (_shape == 0.0 && _invScale == 0.0) if (_shape == 0.0 && _invScale == 0.0)
{ {
return 0.0; return 0.0;
} }
return SpecialFunctions.GammaLowerRegularized(_shape, x * _invScale); return SpecialFunctions.GammaLowerRegularized(_shape, x * _invScale);
} }
@ -361,12 +343,12 @@ namespace MathNet.Numerics.Distributions
{ {
return _shape; return _shape;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return (_shape - 1.0) / _invScale; return (_shape - 1.0) / _invScale;
} }
} }
@ -376,10 +358,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -387,10 +366,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -398,10 +374,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return double.PositiveInfinity; }
{
return double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -420,12 +393,12 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (_shape == 1.0) if (_shape == 1.0)
{ {
return _invScale * Math.Exp(-_invScale * x); return _invScale * Math.Exp(-_invScale * x);
} }
return Math.Pow(_invScale, _shape) * Math.Pow(x, _shape - 1.0) * Math.Exp(-_invScale * x) / SpecialFunctions.Gamma(_shape); return Math.Pow(_invScale, _shape) * Math.Pow(x, _shape - 1.0) * Math.Exp(-_invScale * x) / SpecialFunctions.Gamma(_shape);
} }
@ -440,39 +413,18 @@ namespace MathNet.Numerics.Distributions
{ {
return x == _shape ? Double.PositiveInfinity : Double.NegativeInfinity; return x == _shape ? Double.PositiveInfinity : Double.NegativeInfinity;
} }
if (_shape == 0.0 && _invScale == 0.0) if (_shape == 0.0 && _invScale == 0.0)
{ {
return Double.NegativeInfinity; return Double.NegativeInfinity;
} }
if (_shape == 1.0) if (_shape == 1.0)
{ {
return Math.Log(_invScale) - (_invScale * x); return Math.Log(_invScale) - (_invScale * x);
} }
return (_shape * Math.Log(_invScale)) + ((_shape - 1.0) * Math.Log(x)) - (_invScale * x) - SpecialFunctions.GammaLn(_shape);
}
/// <summary>
/// Generates a sample from the Erlang distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return DoSample(RandomSource, _shape, _invScale);
}
/// <summary> return (_shape * Math.Log(_invScale)) + ((_shape - 1.0) * Math.Log(x)) - (_invScale * x) - SpecialFunctions.GammaLn(_shape);
/// Generates a sequence of samples from the Erlang distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
while (true)
{
yield return DoSample(RandomSource, _shape, _invScale);
}
} }
#endregion #endregion
@ -487,13 +439,13 @@ namespace MathNet.Numerics.Distributions
/// <param name="shape">The shape of the Gamma distribution.</param> /// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param> /// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <returns>A sample from a Erlang distributed random variable.</returns> /// <returns>A sample from a Erlang distributed random variable.</returns>
private static double DoSample(Random rnd, double shape, double invScale) internal static double SampleUnchecked(Random rnd, double shape, double invScale)
{ {
if (Double.IsPositiveInfinity(invScale)) if (Double.IsPositiveInfinity(invScale))
{ {
return shape; return shape;
} }
var a = shape; var a = shape;
var alphafix = 1.0; var alphafix = 1.0;
@ -530,5 +482,63 @@ namespace MathNet.Numerics.Distributions
} }
} }
} }
/// <summary>
/// Generates a sample from the Erlang distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return SampleUnchecked(RandomSource, _shape, _invScale);
}
/// <summary>
/// Generates a sequence of samples from the Erlang distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
while (true)
{
yield return SampleUnchecked(RandomSource, _shape, _invScale);
}
}
/// <summary>
/// Generates a sample from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rnd, double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, shape, invScale);
}
/// <summary>
/// Generates a sequence of samples from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rnd, double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleUnchecked(rnd, shape, invScale);
}
}
} }
} }

115
src/Numerics/Distributions/Continuous/Exponential.cs

@ -44,12 +44,12 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The lambda parameter of the Exponential distribution. /// The lambda parameter of the Exponential distribution.
/// </summary> /// </summary>
private double _lambda; double _lambda;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="Exponential"/> class. /// Initializes a new instance of the <see cref="Exponential"/> class.
@ -68,7 +68,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="lambda">Lambda parameter.</param> /// <param name="lambda">Lambda parameter.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double lambda) void SetParameters(double lambda)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(lambda)) if (Control.CheckDistributionParameters && !IsValidParameterSet(lambda))
{ {
@ -83,7 +83,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="lambda">Lambda parameter.</param> /// <param name="lambda">Lambda parameter.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double lambda) static bool IsValidParameterSet(double lambda)
{ {
if (lambda < 0) if (lambda < 0)
{ {
@ -103,15 +103,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Lambda public double Lambda
{ {
get get { return _lambda; }
{
return _lambda;
}
set set { SetParameters(value); }
{
SetParameters(value);
}
} }
/// <summary> /// <summary>
@ -130,10 +124,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -151,10 +142,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get get { return 1.0 / _lambda; }
{
return 1.0 / _lambda;
}
} }
/// <summary> /// <summary>
@ -162,10 +150,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get get { return 1.0 / (_lambda * _lambda); }
{
return 1.0 / (_lambda * _lambda);
}
} }
/// <summary> /// <summary>
@ -173,10 +158,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return 1.0 / _lambda; }
{
return 1.0 / _lambda;
}
} }
/// <summary> /// <summary>
@ -184,10 +166,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return 1.0 - Math.Log(_lambda); }
{
return 1.0 - Math.Log(_lambda);
}
} }
/// <summary> /// <summary>
@ -195,10 +174,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Skewness public double Skewness
{ {
get get { return 2.0; }
{
return 2.0;
}
} }
/// <summary> /// <summary>
@ -225,10 +201,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -236,10 +209,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { return Math.Log(2.0) / _lambda; }
{
return Math.Log(2.0) / _lambda;
}
} }
/// <summary> /// <summary>
@ -247,10 +217,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -258,10 +225,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -289,13 +253,32 @@ namespace MathNet.Numerics.Distributions
return Math.Log(_lambda) - (_lambda * x); return Math.Log(_lambda) - (_lambda * x);
} }
#endregion
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="lambda">The lambda parameter of the Exponential distribution.</param>
/// <returns>a random number from the distribution.</returns>
internal static double SampleUnchecked(Random rnd, double lambda)
{
var r = rnd.NextDouble();
while (r == 0.0)
{
r = rnd.NextDouble();
}
return -Math.Log(r) / lambda;
}
/// <summary> /// <summary>
/// Draws a random sample from the distribution. /// Draws a random sample from the distribution.
/// </summary> /// </summary>
/// <returns>A random number from this distribution.</returns> /// <returns>A random number from this distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, _lambda); return SampleUnchecked(RandomSource, _lambda);
} }
/// <summary> /// <summary>
@ -306,7 +289,7 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource, _lambda); yield return SampleUnchecked(RandomSource, _lambda);
} }
} }
@ -322,8 +305,8 @@ namespace MathNet.Numerics.Distributions
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return DoSample(rnd, lambda); return SampleUnchecked(rnd, lambda);
} }
/// <summary> /// <summary>
@ -341,26 +324,8 @@ namespace MathNet.Numerics.Distributions
while (true) while (true)
{ {
yield return DoSample(rnd, lambda); yield return SampleUnchecked(rnd, lambda);
}
}
#endregion
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="lambda">The lambda parameter of the Exponential distribution.</param>
/// <returns>a random number from the distribution.</returns>
private static double DoSample(Random rnd, double lambda)
{
var r = rnd.NextDouble();
while (r == 0.0)
{
r = rnd.NextDouble();
} }
return -Math.Log(r) / lambda;
} }
} }
} }

111
src/Numerics/Distributions/Continuous/FisherSnedecor.cs

@ -44,17 +44,17 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The first parameter - degree of freedom. /// The first parameter - degree of freedom.
/// </summary> /// </summary>
private double _d1; double _d1;
/// <summary> /// <summary>
/// The second parameter - degree of freedom. /// The second parameter - degree of freedom.
/// </summary> /// </summary>
private double _d2; double _d2;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="FisherSnedecor"/> class. /// Initializes a new instance of the <see cref="FisherSnedecor"/> class.
@ -76,7 +76,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="d1">The first parameter - degree of freedom.</param> /// <param name="d1">The first parameter - degree of freedom.</param>
/// <param name="d2">The second parameter - degree of freedom.</param> /// <param name="d2">The second parameter - degree of freedom.</param>
private void SetParameters(double d1, double d2) void SetParameters(double d1, double d2)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(d1, d2)) if (Control.CheckDistributionParameters && !IsValidParameterSet(d1, d2))
{ {
@ -93,7 +93,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="d1">The first parameter - degree of freedom.</param> /// <param name="d1">The first parameter - degree of freedom.</param>
/// <param name="d2">The second parameter - degree of freedom.</param> /// <param name="d2">The second parameter - degree of freedom.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double d1, double d2) static bool IsValidParameterSet(double d1, double d2)
{ {
if (d1 <= 0 || d2 <= 0) if (d1 <= 0 || d2 <= 0)
{ {
@ -113,15 +113,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double DegreeOfFreedom1 public double DegreeOfFreedom1
{ {
get get { return _d1; }
{
return _d1;
}
set set { SetParameters(value, _d2); }
{
SetParameters(value, _d2);
}
} }
/// <summary> /// <summary>
@ -129,15 +123,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double DegreeOfFreedom2 public double DegreeOfFreedom2
{ {
get get { return _d2; }
{
return _d2;
}
set set { SetParameters(_d1, value); }
{
SetParameters(_d1, value);
}
} }
/// <summary> /// <summary>
@ -156,10 +144,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -209,10 +194,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return Math.Sqrt(Variance); }
{
return Math.Sqrt(Variance);
}
} }
/// <summary> /// <summary>
@ -220,10 +202,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -277,10 +256,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -288,10 +264,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -299,10 +272,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return double.PositiveInfinity; }
{
return double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -325,13 +295,27 @@ namespace MathNet.Numerics.Distributions
return Math.Log(Density(x)); return Math.Log(Density(x));
} }
#endregion
/// <summary>
/// Generates one sample from the <c>FisherSnedecor</c> distribution without parameter checking.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="d1">The first parameter - degree of freedom.</param>
/// <param name="d2">The second parameter - degree of freedom.</param>
/// <returns>a <c>FisherSnedecor</c> distributed random number.</returns>
internal static double SampleUnchecked(Random rnd, double d1, double d2)
{
return (ChiSquare.Sample(rnd, d1) / d1) / (ChiSquare.Sample(rnd, d2) / d2);
}
/// <summary> /// <summary>
/// Generates a sample from the <c>FisherSnedecor</c> distribution. /// Generates a sample from the <c>FisherSnedecor</c> distribution.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, _d1, _d2); return SampleUnchecked(RandomSource, _d1, _d2);
} }
/// <summary> /// <summary>
@ -342,22 +326,45 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource, _d1, _d2); yield return SampleUnchecked(RandomSource, _d1, _d2);
} }
} }
#endregion /// <summary>
/// Generates a sample from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="d1">The first parameter - degree of freedom.</param>
/// <param name="d2">The second parameter - degree of freedom.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rnd, double d1, double d2)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(d1, d2))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, d1, d2);
}
/// <summary> /// <summary>
/// Generates one sample from the <c>FisherSnedecor</c> distribution without parameter checking. /// Generates a sequence of samples from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="d1">The first parameter - degree of freedom.</param> /// <param name="d1">The first parameter - degree of freedom.</param>
/// <param name="d2">The second parameter - degree of freedom.</param> /// <param name="d2">The second parameter - degree of freedom.</param>
/// <returns>a <c>FisherSnedecor</c> distributed random number.</returns> /// <returns>a sequence of samples from the distribution.</returns>
private static double DoSample(Random rnd, double d1, double d2) public static IEnumerable<double> Samples(Random rnd, double d1, double d2)
{ {
return (ChiSquare.Sample(rnd, d1) / d1) / (ChiSquare.Sample(rnd, d2) / d2); if (Control.CheckDistributionParameters && !IsValidParameterSet(d1, d2))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleUnchecked(rnd, d1, d2);
}
} }
} }
} }

213
src/Numerics/Distributions/Continuous/Gamma.cs

@ -52,17 +52,17 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Gamma shape parameter. /// Gamma shape parameter.
/// </summary> /// </summary>
private double _shape; double _shape;
/// <summary> /// <summary>
/// Gamma inverse scale parameter. /// Gamma inverse scale parameter.
/// </summary> /// </summary>
private double _invScale; double _invScale;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the Gamma class. /// Initializes a new instance of the Gamma class.
@ -114,7 +114,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="shape">The shape of the Gamma distribution.</param> /// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param> /// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double shape, double invScale) static bool IsValidParameterSet(double shape, double invScale)
{ {
if (shape < 0.0 || invScale < 0.0 || Double.IsNaN(shape) || Double.IsNaN(invScale)) if (shape < 0.0 || invScale < 0.0 || Double.IsNaN(shape) || Double.IsNaN(invScale))
{ {
@ -130,7 +130,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="shape">The shape of the Gamma distribution.</param> /// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param> /// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double shape, double invScale) void SetParameters(double shape, double invScale)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale)) if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{ {
@ -146,15 +146,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Shape public double Shape
{ {
get get { return _shape; }
{
return _shape;
}
set set { SetParameters(value, _invScale); }
{
SetParameters(value, _invScale);
}
} }
/// <summary> /// <summary>
@ -162,10 +156,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return 1.0 / _invScale; }
{
return 1.0 / _invScale;
}
set set
{ {
@ -185,15 +176,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double InvScale public double InvScale
{ {
get get { return _invScale; }
{
return _invScale;
}
set set { SetParameters(_shape, value); }
{
SetParameters(_shape, value);
}
} }
#region IDistribution implementation #region IDistribution implementation
@ -203,10 +188,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -230,12 +212,12 @@ namespace MathNet.Numerics.Distributions
{ {
return _shape; return _shape;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return _shape / _invScale; return _shape / _invScale;
} }
} }
@ -251,12 +233,12 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return _shape / (_invScale * _invScale); return _shape / (_invScale * _invScale);
} }
} }
@ -272,12 +254,12 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return Math.Sqrt(_shape / (_invScale * _invScale)); return Math.Sqrt(_shape / (_invScale * _invScale));
} }
} }
@ -293,12 +275,12 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return _shape - Math.Log(_invScale) + SpecialFunctions.GammaLn(_shape) + ((1.0 - _shape) * SpecialFunctions.DiGamma(_shape)); return _shape - Math.Log(_invScale) + SpecialFunctions.GammaLn(_shape) + ((1.0 - _shape) * SpecialFunctions.DiGamma(_shape));
} }
} }
@ -314,12 +296,12 @@ namespace MathNet.Numerics.Distributions
{ {
return 0.0; return 0.0;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return 2.0 / Math.Sqrt(_shape); return 2.0 / Math.Sqrt(_shape);
} }
} }
@ -339,12 +321,12 @@ namespace MathNet.Numerics.Distributions
{ {
return _shape; return _shape;
} }
if (_invScale == 0.0 && _shape == 0.0) if (_invScale == 0.0 && _shape == 0.0)
{ {
return Double.NaN; return Double.NaN;
} }
return (_shape - 1.0) / _invScale; return (_shape - 1.0) / _invScale;
} }
} }
@ -354,10 +336,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -365,10 +344,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -376,10 +352,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -393,17 +366,17 @@ namespace MathNet.Numerics.Distributions
{ {
return x == _shape ? Double.PositiveInfinity : 0.0; return x == _shape ? Double.PositiveInfinity : 0.0;
} }
if (_shape == 0.0 && _invScale == 0.0) if (_shape == 0.0 && _invScale == 0.0)
{ {
return 0.0; return 0.0;
} }
if (_shape == 1.0) if (_shape == 1.0)
{ {
return _invScale * Math.Exp(-_invScale * x); return _invScale * Math.Exp(-_invScale * x);
} }
return Math.Pow(_invScale, _shape) * Math.Pow(x, _shape - 1.0) * Math.Exp(-_invScale * x) / SpecialFunctions.Gamma(_shape); return Math.Pow(_invScale, _shape) * Math.Pow(x, _shape - 1.0) * Math.Exp(-_invScale * x) / SpecialFunctions.Gamma(_shape);
} }
@ -418,17 +391,17 @@ namespace MathNet.Numerics.Distributions
{ {
return x == _shape ? Double.PositiveInfinity : Double.NegativeInfinity; return x == _shape ? Double.PositiveInfinity : Double.NegativeInfinity;
} }
if (_shape == 0.0 && _invScale == 0.0) if (_shape == 0.0 && _invScale == 0.0)
{ {
return Double.NegativeInfinity; return Double.NegativeInfinity;
} }
if (_shape == 1.0) if (_shape == 1.0)
{ {
return Math.Log(_invScale) - (_invScale * x); return Math.Log(_invScale) - (_invScale * x);
} }
return (_shape * Math.Log(_invScale)) + ((_shape - 1.0) * Math.Log(x)) - (_invScale * x) - SpecialFunctions.GammaLn(_shape); return (_shape * Math.Log(_invScale)) + ((_shape - 1.0) * Math.Log(x)) - (_invScale * x) - SpecialFunctions.GammaLn(_shape);
} }
@ -443,75 +416,17 @@ namespace MathNet.Numerics.Distributions
{ {
return x >= _shape ? 1.0 : 0.0; return x >= _shape ? 1.0 : 0.0;
} }
if (_shape == 0.0 && _invScale == 0.0) if (_shape == 0.0 && _invScale == 0.0)
{ {
return 0.0; return 0.0;
} }
return SpecialFunctions.GammaLowerRegularized(_shape, x * _invScale);
}
/// <summary> return SpecialFunctions.GammaLowerRegularized(_shape, x * _invScale);
/// Generates a sample from the Gamma distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return SampleGamma(RandomSource, _shape, _invScale);
}
/// <summary>
/// Generates a sequence of samples from the Gamma distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
while (true)
{
yield return SampleGamma(RandomSource, _shape, _invScale);
}
} }
#endregion #endregion
/// <summary>
/// Generates a sample from the Gamma distribution.
/// </summary>
/// <param name="rng">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution from which to generate samples.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution from which to generate samples.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rng, double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleGamma(rng, shape, invScale);
}
/// <summary>
/// Generates a sequence of samples from the Gamma distribution.
/// </summary>
/// <param name="rng">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution from which to generate samples.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution from which to generate samples.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rng, double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleGamma(rng, shape, invScale);
}
}
/// <summary> /// <summary>
/// <para>Sampling implementation based on: /// <para>Sampling implementation based on:
/// "A Simple Method for Generating Gamma Variables" - Marsaglia &amp; Tsang /// "A Simple Method for Generating Gamma Variables" - Marsaglia &amp; Tsang
@ -522,7 +437,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="shape">The shape of the Gamma distribution.</param> /// <param name="shape">The shape of the Gamma distribution.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution.</param> /// <param name="invScale">The inverse scale of the Gamma distribution.</param>
/// <returns>A sample from a Gamma distributed random variable.</returns> /// <returns>A sample from a Gamma distributed random variable.</returns>
internal static double SampleGamma(Random rnd, double shape, double invScale) internal static double SampleUnchecked(Random rnd, double shape, double invScale)
{ {
if (Double.IsPositiveInfinity(invScale)) if (Double.IsPositiveInfinity(invScale))
{ {
@ -565,5 +480,63 @@ namespace MathNet.Numerics.Distributions
} }
} }
} }
/// <summary>
/// Generates a sample from the Gamma distribution.
/// </summary>
/// <returns>a sample from the distribution.</returns>
public double Sample()
{
return SampleUnchecked(RandomSource, _shape, _invScale);
}
/// <summary>
/// Generates a sequence of samples from the Gamma distribution.
/// </summary>
/// <returns>a sequence of samples from the distribution.</returns>
public IEnumerable<double> Samples()
{
while (true)
{
yield return SampleUnchecked(RandomSource, _shape, _invScale);
}
}
/// <summary>
/// Generates a sample from the Gamma distribution.
/// </summary>
/// <param name="rng">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution from which to generate samples.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution from which to generate samples.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rng, double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rng, shape, invScale);
}
/// <summary>
/// Generates a sequence of samples from the Gamma distribution.
/// </summary>
/// <param name="rng">The random number generator to use.</param>
/// <param name="shape">The shape of the Gamma distribution from which to generate samples.</param>
/// <param name="invScale">The inverse scale of the Gamma distribution from which to generate samples.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rng, double shape, double invScale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, invScale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleUnchecked(rng, shape, invScale);
}
}
} }
} }

126
src/Numerics/Distributions/Continuous/InverseGamma.cs

@ -45,17 +45,17 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Inverse Gamma shape parameter. /// Inverse Gamma shape parameter.
/// </summary> /// </summary>
private double _shape; double _shape;
/// <summary> /// <summary>
/// Inverse Gamma scale parameter scale. /// Inverse Gamma scale parameter scale.
/// </summary> /// </summary>
private double _scale; double _scale;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="InverseGamma"/> class. /// Initializes a new instance of the <see cref="InverseGamma"/> class.
@ -82,7 +82,7 @@ namespace MathNet.Numerics.Distributions
/// The scale (beta) parameter of the inverse Gamma distribution. /// The scale (beta) parameter of the inverse Gamma distribution.
/// </param> /// </param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double shape, double scale) void SetParameters(double shape, double scale)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, scale)) if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, scale))
{ {
@ -103,7 +103,7 @@ namespace MathNet.Numerics.Distributions
/// The scale (beta) parameter of the inverse Gamma distribution. /// The scale (beta) parameter of the inverse Gamma distribution.
/// </param> /// </param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double shape, double scale) static bool IsValidParameterSet(double shape, double scale)
{ {
if (shape <= 0 || scale <= 0) if (shape <= 0 || scale <= 0)
{ {
@ -123,15 +123,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Shape public double Shape
{ {
get get { return _shape; }
{
return _shape;
}
set set { SetParameters(value, _scale); }
{
SetParameters(value, _scale);
}
} }
/// <summary> /// <summary>
@ -139,15 +133,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return _scale; }
{
return _scale;
}
set set { SetParameters(_shape, value); }
{
SetParameters(_shape, value);
}
} }
/// <summary> /// <summary>
@ -166,10 +154,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -219,10 +204,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return _scale / (Math.Abs(_shape - 1.0) * Math.Sqrt(_shape - 2.0)); }
{
return _scale / (Math.Abs(_shape - 1.0) * Math.Sqrt(_shape - 2.0));
}
} }
/// <summary> /// <summary>
@ -230,10 +212,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return _shape + Math.Log(_scale) + SpecialFunctions.GammaLn(_shape) - ((1 + _shape) * SpecialFunctions.DiGamma(_shape)); }
{
return _shape + Math.Log(_scale) + SpecialFunctions.GammaLn(_shape) - ((1 + _shape) * SpecialFunctions.DiGamma(_shape));
}
} }
/// <summary> /// <summary>
@ -271,10 +250,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return _scale / (_shape + 1.0); }
{
return _scale / (_shape + 1.0);
}
} }
/// <summary> /// <summary>
@ -283,10 +259,7 @@ namespace MathNet.Numerics.Distributions
/// <remarks>Throws <see cref="NotSupportedException"/>.</remarks> /// <remarks>Throws <see cref="NotSupportedException"/>.</remarks>
public double Median public double Median
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -294,10 +267,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -305,10 +275,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -336,13 +303,27 @@ namespace MathNet.Numerics.Distributions
return Math.Log(Density(x)); return Math.Log(Density(x));
} }
#endregion
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="shape">The shape (alpha) parameter of the inverse Gamma distribution.</param>
/// <param name="scale">The scale (beta) parameter of the inverse Gamma distribution.</param>
/// <returns>a random number from the distribution.</returns>
internal static double SampleUnchecked(Random rnd, double shape, double scale)
{
return 1.0 / Gamma.Sample(rnd, shape, scale);
}
/// <summary> /// <summary>
/// Draws a random sample from the distribution. /// Draws a random sample from the distribution.
/// </summary> /// </summary>
/// <returns>A random number from this distribution.</returns> /// <returns>A random number from this distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, _shape, _scale); return SampleUnchecked(RandomSource, _shape, _scale);
} }
/// <summary> /// <summary>
@ -353,26 +334,45 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource, _shape, _scale); yield return SampleUnchecked(RandomSource, _shape, _scale);
} }
} }
#endregion
/// <summary> /// <summary>
/// Samples the distribution. /// Generates a sample from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="shape"> /// <param name="shape">The shape (alpha) parameter of the inverse Gamma distribution.</param>
/// The shape (alpha) parameter of the inverse Gamma distribution. /// <param name="scale">The scale (beta) parameter of the inverse Gamma distribution.</param>
/// </param> /// <returns>a sample from the distribution.</returns>
/// <param name="scale"> public static double Sample(Random rnd, double shape, double scale)
/// The scale (beta) parameter of the inverse Gamma distribution.
/// </param>
/// <returns>a random number from the distribution.</returns>
private static double DoSample(Random rnd, double shape, double scale)
{ {
return 1.0 / Gamma.Sample(rnd, shape, scale); if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, scale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, shape, scale);
}
/// <summary>
/// Generates a sequence of samples from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="shape">The shape (alpha) parameter of the inverse Gamma distribution.</param>
/// <param name="scale">The scale (beta) parameter of the inverse Gamma distribution.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rnd, double shape, double scale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, scale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleUnchecked(rnd, shape, scale);
}
} }
} }
} }

135
src/Numerics/Distributions/Continuous/Laplace.cs

@ -46,27 +46,21 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The scale of the Laplace distribution. /// The scale of the Laplace distribution.
/// </summary> /// </summary>
private double _scale; double _scale;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Gets or sets the location of the Laplace distribution. /// Gets or sets the location of the Laplace distribution.
/// </summary> /// </summary>
public double Location public double Location
{ {
get get { return Mean; }
{
return Mean;
}
set set { SetParameters(value, _scale); }
{
SetParameters(value, _scale);
}
} }
/// <summary> /// <summary>
@ -74,21 +68,16 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return _scale; }
{
return _scale;
}
set set { SetParameters(Mean, value); }
{
SetParameters(Mean, value);
}
} }
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="Laplace"/> class (location = 0, scale = 1). /// Initializes a new instance of the <see cref="Laplace"/> class (location = 0, scale = 1).
/// </summary> /// </summary>
public Laplace() : this(0.0, 1.0) public Laplace()
: this(0.0, 1.0)
{ {
} }
@ -116,7 +105,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="location">The location for the Laplace distribution.</param> /// <param name="location">The location for the Laplace distribution.</param>
/// <param name="scale">The scale for the Laplace distribution.</param> /// <param name="scale">The scale for the Laplace distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double location, double scale) void SetParameters(double location, double scale)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale)) if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale))
{ {
@ -133,7 +122,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="location">The location for the Laplace distribution.</param> /// <param name="location">The location for the Laplace distribution.</param>
/// <param name="scale">The scale for the Laplace distribution.</param> /// <param name="scale">The scale for the Laplace distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double location, double scale) static bool IsValidParameterSet(double location, double scale)
{ {
if (scale <= 0) if (scale <= 0)
{ {
@ -164,10 +153,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -183,21 +169,14 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Gets the mean of the distribution. /// Gets the mean of the distribution.
/// </summary> /// </summary>
public double Mean public double Mean { get; private set; }
{
get;
private set;
}
/// <summary> /// <summary>
/// Gets the variance of the distribution. /// Gets the variance of the distribution.
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get get { return 2.0 * _scale * _scale; }
{
return 2.0 * _scale * _scale;
}
} }
/// <summary> /// <summary>
@ -205,10 +184,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return Math.Sqrt(2.0) * _scale; }
{
return Math.Sqrt(2.0) * _scale;
}
} }
/// <summary> /// <summary>
@ -216,10 +192,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return Math.Log(2.0 * Constants.E * _scale); }
{
return Math.Log(2.0 * Constants.E * _scale);
}
} }
/// <summary> /// <summary>
@ -227,10 +200,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Skewness public double Skewness
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -252,10 +222,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return Mean; }
{
return Mean;
}
} }
/// <summary> /// <summary>
@ -263,10 +230,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { return Mean; }
{
return Mean;
}
} }
/// <summary> /// <summary>
@ -274,10 +238,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return Double.NegativeInfinity; }
{
return Double.NegativeInfinity;
}
} }
/// <summary> /// <summary>
@ -285,10 +246,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -311,13 +269,28 @@ namespace MathNet.Numerics.Distributions
return Math.Log(Density(x)); return Math.Log(Density(x));
} }
#endregion
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="location">The location shape parameter.</param>
/// <param name="scale">The scale parameter.</param>
/// <returns>a random number from the distribution.</returns>
internal static double SampleUnchecked(Random rnd, double location, double scale)
{
var u = rnd.NextDouble() - 0.5;
return location - (scale * Math.Sign(u) * Math.Log(1.0 - (2.0 * Math.Abs(u))));
}
/// <summary> /// <summary>
/// Samples a Laplace distributed random variable. /// Samples a Laplace distributed random variable.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, Mean, _scale); return SampleUnchecked(RandomSource, Mean, _scale);
} }
/// <summary> /// <summary>
@ -328,23 +301,45 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource, Mean, _scale); yield return SampleUnchecked(RandomSource, Mean, _scale);
} }
} }
#endregion /// <summary>
/// Generates a sample from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="location">The location shape parameter.</param>
/// <param name="scale">The scale parameter.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rnd, double location, double scale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, location, scale);
}
/// <summary> /// <summary>
/// Samples the distribution. /// Generates a sequence of samples from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="location">The location shape parameter.</param> /// <param name="location">The location shape parameter.</param>
/// <param name="scale">The scale parameter.</param> /// <param name="scale">The scale parameter.</param>
/// <returns>a random number from the distribution.</returns> /// <returns>a sequence of samples from the distribution.</returns>
private static double DoSample(Random rnd, double location, double scale) public static IEnumerable<double> Samples(Random rnd, double location, double scale)
{ {
var u = rnd.NextDouble() - 0.5; if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale))
return location - (scale * Math.Sign(u) * Math.Log(1.0 - (2.0 * Math.Abs(u)))); {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleUnchecked(rnd, location, scale);
}
} }
} }
} }

77
src/Numerics/Distributions/Continuous/LogNormal.cs

@ -44,17 +44,17 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Keeps track of the mu of the logarithm of the log-log-normal distribution. /// Keeps track of the mu of the logarithm of the log-log-normal distribution.
/// </summary> /// </summary>
private double _mu; double _mu;
/// <summary> /// <summary>
/// Keeps track of the standard deviation of the logarithm of the log-log-normal distribution. /// Keeps track of the standard deviation of the logarithm of the log-log-normal distribution.
/// </summary> /// </summary>
private double _sigma; double _sigma;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="LogNormal"/> class. /// Initializes a new instance of the <see cref="LogNormal"/> class.
@ -88,7 +88,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="mu">The mu of the logarithm of the distribution.</param> /// <param name="mu">The mu of the logarithm of the distribution.</param>
/// <param name="sigma">The standard deviation of the logarithm of the distribution.</param> /// <param name="sigma">The standard deviation of the logarithm of the distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double mu, double sigma) static bool IsValidParameterSet(double mu, double sigma)
{ {
if (sigma < 0.0 || Double.IsNaN(mu) || Double.IsNaN(mu) || Double.IsNaN(sigma)) if (sigma < 0.0 || Double.IsNaN(mu) || Double.IsNaN(mu) || Double.IsNaN(sigma))
{ {
@ -104,7 +104,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="mu">The mu of the logarithm of the distribution.</param> /// <param name="mu">The mu of the logarithm of the distribution.</param>
/// <param name="sigma">The standard deviation of the logarithm of the distribution.</param> /// <param name="sigma">The standard deviation of the logarithm of the distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double mu, double sigma) void SetParameters(double mu, double sigma)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(mu, sigma)) if (Control.CheckDistributionParameters && !IsValidParameterSet(mu, sigma))
{ {
@ -120,15 +120,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mu public double Mu
{ {
get get { return _mu; }
{
return _mu;
}
set set { SetParameters(value, _sigma); }
{
SetParameters(value, _sigma);
}
} }
/// <summary> /// <summary>
@ -136,15 +130,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Sigma public double Sigma
{ {
get get { return _sigma; }
{
return _sigma;
}
set set { SetParameters(_mu, value); }
{
SetParameters(_mu, value);
}
} }
#region IDistribution implementation #region IDistribution implementation
@ -154,10 +142,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -175,10 +160,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get get { return Math.Exp(_mu + (_sigma * _sigma / 2.0)); }
{
return Math.Exp(_mu + (_sigma * _sigma / 2.0));
}
} }
/// <summary> /// <summary>
@ -210,10 +192,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return 0.5 + Math.Log(_sigma) + _mu + Constants.LogSqrt2Pi; }
{
return 0.5 + Math.Log(_sigma) + _mu + Constants.LogSqrt2Pi;
}
} }
/// <summary> /// <summary>
@ -237,10 +216,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return Math.Exp(_mu - (_sigma * _sigma)); }
{
return Math.Exp(_mu - (_sigma * _sigma));
}
} }
/// <summary> /// <summary>
@ -248,10 +224,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { return Math.Exp(_mu); }
{
return Math.Exp(_mu);
}
} }
/// <summary> /// <summary>
@ -259,10 +232,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -270,10 +240,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -323,13 +290,15 @@ namespace MathNet.Numerics.Distributions
return 0.5 * (1.0 + SpecialFunctions.Erf((Math.Log(x) - _mu) / (_sigma * Constants.Sqrt2))); return 0.5 * (1.0 + SpecialFunctions.Erf((Math.Log(x) - _mu) / (_sigma * Constants.Sqrt2)));
} }
#endregion
/// <summary> /// <summary>
/// Generates a sample from the log-normal distribution using the <i>Box-Muller</i> algorithm. /// Generates a sample from the log-normal distribution using the <i>Box-Muller</i> algorithm.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return Math.Exp(_mu + (_sigma * Normal.SampleBoxMuller(RandomSource).Item1)); return Math.Exp(Normal.SampleUnchecked(RandomSource, _mu, _sigma));
} }
/// <summary> /// <summary>
@ -340,14 +309,12 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
var sample = Normal.SampleBoxMuller(RandomSource); var sample = Normal.SampleUncheckedBoxMuller(RandomSource);
yield return Math.Exp(_mu + (_sigma * sample.Item1)); yield return Math.Exp(_mu + (_sigma * sample.Item1));
yield return Math.Exp(_mu + (_sigma * sample.Item2)); yield return Math.Exp(_mu + (_sigma * sample.Item2));
} }
} }
#endregion
/// <summary> /// <summary>
/// Generates a sample from the log-normal distribution using the <i>Box-Muller</i> algorithm. /// Generates a sample from the log-normal distribution using the <i>Box-Muller</i> algorithm.
/// </summary> /// </summary>
@ -362,7 +329,7 @@ namespace MathNet.Numerics.Distributions
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return Math.Exp(mu + (sigma * Normal.SampleBoxMuller(rng).Item1)); return Math.Exp(Normal.SampleUnchecked(rng, mu, sigma));
} }
/// <summary> /// <summary>
@ -381,7 +348,7 @@ namespace MathNet.Numerics.Distributions
while (true) while (true)
{ {
var sample = Normal.SampleBoxMuller(rng); var sample = Normal.SampleUncheckedBoxMuller(rng);
yield return Math.Exp(mu + (sigma * sample.Item1)); yield return Math.Exp(mu + (sigma * sample.Item1));
yield return Math.Exp(mu + (sigma * sample.Item2)); yield return Math.Exp(mu + (sigma * sample.Item2));
} }

179
src/Numerics/Distributions/Continuous/Normal.cs

@ -44,24 +44,25 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Keeps track of the mean of the normal distribution. /// Keeps track of the mean of the normal distribution.
/// </summary> /// </summary>
private double _mean; double _mean;
/// <summary> /// <summary>
/// Keeps track of the standard deviation of the normal distribution. /// Keeps track of the standard deviation of the normal distribution.
/// </summary> /// </summary>
private double _stdDev; double _stdDev;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the Normal class. This is a normal distribution with mean 0.0 /// Initializes a new instance of the Normal class. This is a normal distribution with mean 0.0
/// and standard deviation 1.0. The distribution will /// and standard deviation 1.0. The distribution will
/// be initialized with the default <seealso cref="System.Random"/> random number generator. /// be initialized with the default <seealso cref="System.Random"/> random number generator.
/// </summary> /// </summary>
public Normal() : this(0.0, 1.0) public Normal()
: this(0.0, 1.0)
{ {
} }
@ -128,7 +129,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="mean">The mean of the normal distribution.</param> /// <param name="mean">The mean of the normal distribution.</param>
/// <param name="stddev">The standard deviation of the normal distribution.</param> /// <param name="stddev">The standard deviation of the normal distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double mean, double stddev) static bool IsValidParameterSet(double mean, double stddev)
{ {
if (stddev < 0.0 || Double.IsNaN(mean) || Double.IsNaN(stddev)) if (stddev < 0.0 || Double.IsNaN(mean) || Double.IsNaN(stddev))
{ {
@ -144,7 +145,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="mean">The mean of the normal distribution.</param> /// <param name="mean">The mean of the normal distribution.</param>
/// <param name="stddev">The standard deviation of the normal distribution.</param> /// <param name="stddev">The standard deviation of the normal distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double mean, double stddev) void SetParameters(double mean, double stddev)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev)) if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev))
{ {
@ -160,10 +161,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Precision public double Precision
{ {
get get { return 1.0 / (_stdDev * _stdDev); }
{
return 1.0 / (_stdDev * _stdDev);
}
set set
{ {
@ -186,10 +184,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -207,15 +202,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get get { return _mean; }
{
return _mean;
}
set set { SetParameters(value, _stdDev); }
{
SetParameters(value, _stdDev);
}
} }
/// <summary> /// <summary>
@ -223,15 +212,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get get { return _stdDev * _stdDev; }
{
return _stdDev * _stdDev;
}
set set { SetParameters(_mean, Math.Sqrt(value)); }
{
SetParameters(_mean, Math.Sqrt(value));
}
} }
/// <summary> /// <summary>
@ -239,15 +222,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return _stdDev; }
{
return _stdDev;
}
set set { SetParameters(_mean, value); }
{
SetParameters(_mean, value);
}
} }
/// <summary> /// <summary>
@ -255,10 +232,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return Math.Log(_stdDev) + Constants.LogSqrt2PiE; }
{
return Math.Log(_stdDev) + Constants.LogSqrt2PiE;
}
} }
/// <summary> /// <summary>
@ -266,10 +240,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Skewness public double Skewness
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
#endregion #endregion
@ -281,10 +252,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return _mean; }
{
return _mean;
}
} }
/// <summary> /// <summary>
@ -292,10 +260,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { return _mean; }
{
return _mean;
}
} }
/// <summary> /// <summary>
@ -303,10 +268,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return Double.NegativeInfinity; }
{
return Double.NegativeInfinity;
}
} }
/// <summary> /// <summary>
@ -314,10 +276,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -388,13 +347,60 @@ namespace MathNet.Numerics.Distributions
return CumulativeDistribution(_mean, _stdDev, x); return CumulativeDistribution(_mean, _stdDev, x);
} }
#endregion
/// <summary>
/// Computes the inverse cumulative distribution function of the normal distribution.
/// </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>
public double InverseCumulativeDistribution(double p)
{
return _mean - (_stdDev * Math.Sqrt(2.0) * SpecialFunctions.ErfcInv(2.0 * p));
}
/// <summary>
/// Samples a pair of standard normal distributed random variables using the <i>Box-Muller</i> algorithm.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <returns>a pair of random numbers from the standard normal distribution.</returns>
internal static Tuple<double, double> SampleUncheckedBoxMuller(Random rnd)
{
var v1 = (2.0 * rnd.NextDouble()) - 1.0;
var v2 = (2.0 * rnd.NextDouble()) - 1.0;
var r = (v1 * v1) + (v2 * v2);
while (r >= 1.0 || r == 0.0)
{
v1 = (2.0 * rnd.NextDouble()) - 1.0;
v2 = (2.0 * rnd.NextDouble()) - 1.0;
r = (v1 * v1) + (v2 * v2);
}
var fac = Math.Sqrt(-2.0 * Math.Log(r) / r);
return new Tuple<double, double>(v1 * fac, v2 * fac);
}
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="mean">The mean of the normal distribution from which to generate samples.</param>
/// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param>
/// <returns>a random number from the distribution.</returns>
internal static double SampleUnchecked(Random rnd, double mean, double stddev)
{
return mean + (stddev * SampleUncheckedBoxMuller(rnd).Item1);
}
/// <summary> /// <summary>
/// Generates a sample from the normal distribution using the <i>Box-Muller</i> algorithm. /// Generates a sample from the normal distribution using the <i>Box-Muller</i> algorithm.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return _mean + (_stdDev * SampleBoxMuller(RandomSource).Item1); return SampleUnchecked(RandomSource, _mean, _stdDev);
} }
/// <summary> /// <summary>
@ -405,49 +411,37 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
var sample = SampleBoxMuller(RandomSource); var sample = SampleUncheckedBoxMuller(RandomSource);
yield return _mean + (_stdDev * sample.Item1); yield return _mean + (_stdDev * sample.Item1);
yield return _mean + (_stdDev * sample.Item2); yield return _mean + (_stdDev * sample.Item2);
} }
} }
#endregion
/// <summary>
/// Computes the inverse cumulative distribution function of the normal distribution.
/// </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>
public double InverseCumulativeDistribution(double p)
{
return _mean - (_stdDev * Math.Sqrt(2.0) * SpecialFunctions.ErfcInv(2.0 * p));
}
/// <summary> /// <summary>
/// Generates a sample from the normal distribution using the <i>Box-Muller</i> algorithm. /// Generates a sample from the normal distribution using the <i>Box-Muller</i> algorithm.
/// </summary> /// </summary>
/// <param name="rng">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="mean">The mean of the normal distribution from which to generate samples.</param> /// <param name="mean">The mean of the normal distribution from which to generate samples.</param>
/// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param> /// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public static double Sample(Random rng, double mean, double stddev) public static double Sample(Random rnd, double mean, double stddev)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev)) if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev))
{ {
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return mean + (stddev * SampleBoxMuller(rng).Item1); return SampleUnchecked(rnd, mean, stddev);
} }
/// <summary> /// <summary>
/// Generates a sequence of samples from the normal distribution using the <i>Box-Muller</i> algorithm. /// Generates a sequence of samples from the normal distribution using the <i>Box-Muller</i> algorithm.
/// </summary> /// </summary>
/// <param name="rng">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="mean">The mean of the normal distribution from which to generate samples.</param> /// <param name="mean">The mean of the normal distribution from which to generate samples.</param>
/// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param> /// <param name="stddev">The standard deviation of the normal distribution from which to generate samples.</param>
/// <returns>a sequence of samples from the distribution.</returns> /// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rng, double mean, double stddev) public static IEnumerable<double> Samples(Random rnd, double mean, double stddev)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev)) if (Control.CheckDistributionParameters && !IsValidParameterSet(mean, stddev))
{ {
@ -456,31 +450,10 @@ namespace MathNet.Numerics.Distributions
while (true) while (true)
{ {
var sample = SampleBoxMuller(rng); var sample = SampleUncheckedBoxMuller(rnd);
yield return mean + (stddev * sample.Item1); yield return mean + (stddev * sample.Item1);
yield return mean + (stddev * sample.Item2); yield return mean + (stddev * sample.Item2);
} }
} }
/// <summary>
/// Samples a pair of standard normal distributed random variables using the <i>Box-Muller</i> algorithm.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <returns>a pair of random numbers from the standard normal distribution.</returns>
internal static Tuple<double, double> SampleBoxMuller(Random rnd)
{
var v1 = (2.0 * rnd.NextDouble()) - 1.0;
var v2 = (2.0 * rnd.NextDouble()) - 1.0;
var r = (v1 * v1) + (v2 * v2);
while (r >= 1.0 || r == 0.0)
{
v1 = (2.0 * rnd.NextDouble()) - 1.0;
v2 = (2.0 * rnd.NextDouble()) - 1.0;
r = (v1 * v1) + (v2 * v2);
}
var fac = Math.Sqrt(-2.0 * Math.Log(r) / r);
return new Tuple<double, double>(v1 * fac, v2 * fac);
}
} }
} }

123
src/Numerics/Distributions/Continuous/Pareto.cs

@ -46,17 +46,17 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The scale parameter of the distribution. /// The scale parameter of the distribution.
/// </summary> /// </summary>
private double _scale; double _scale;
/// <summary> /// <summary>
/// The shape parameter of the distribution. /// The shape parameter of the distribution.
/// </summary> /// </summary>
private double _shape; double _shape;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="Pareto"/> class. /// Initializes a new instance of the <see cref="Pareto"/> class.
@ -82,7 +82,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="scale">The scale parameter of the distribution.</param> /// <param name="scale">The scale parameter of the distribution.</param>
/// <param name="shape">The shape parameter of the distribution.</param> /// <param name="shape">The shape parameter of the distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double scale, double shape) void SetParameters(double scale, double shape)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(scale, shape)) if (Control.CheckDistributionParameters && !IsValidParameterSet(scale, shape))
{ {
@ -99,7 +99,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="scale">The scale parameter of the distribution.</param> /// <param name="scale">The scale parameter of the distribution.</param>
/// <param name="shape">The shape parameter of the distribution.</param> /// <param name="shape">The shape parameter of the distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double scale, double shape) static bool IsValidParameterSet(double scale, double shape)
{ {
if (scale <= 0 || shape <= 0) if (scale <= 0 || shape <= 0)
{ {
@ -119,15 +119,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return _scale; }
{
return _scale;
}
set set { SetParameters(value, _shape); }
{
SetParameters(value, _shape);
}
} }
/// <summary> /// <summary>
@ -135,15 +129,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Shape public double Shape
{ {
get get { return _shape; }
{
return _shape;
}
set set { SetParameters(_scale, value); }
{
SetParameters(_scale, value);
}
} }
/// <summary> /// <summary>
@ -162,10 +150,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -215,10 +200,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return (_scale * Math.Sqrt(_shape)) / (Math.Abs(_shape - 1.0) * Math.Sqrt(_shape - 2.0)); }
{
return (_scale * Math.Sqrt(_shape)) / (Math.Abs(_shape - 1.0) * Math.Sqrt(_shape - 2.0));
}
} }
/// <summary> /// <summary>
@ -226,10 +208,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return Math.Log(_shape / _scale) - (1.0 / _shape) - 1.0; }
{
return Math.Log(_shape / _scale) - (1.0 / _shape) - 1.0;
}
} }
/// <summary> /// <summary>
@ -237,10 +216,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Skewness public double Skewness
{ {
get get { return (2.0 * (_shape + 1.0) / (_shape - 3.0)) * Math.Sqrt((_shape - 2.0) / _shape); }
{
return (2.0 * (_shape + 1.0) / (_shape - 3.0)) * Math.Sqrt((_shape - 2.0) / _shape);
}
} }
/// <summary> /// <summary>
@ -262,10 +238,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return _scale; }
{
return _scale;
}
} }
/// <summary> /// <summary>
@ -273,10 +246,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { return _scale * Math.Pow(2.0, 1.0 / _shape); }
{
return _scale * Math.Pow(2.0, 1.0 / _shape);
}
} }
/// <summary> /// <summary>
@ -284,10 +254,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return _scale; }
{
return _scale;
}
} }
/// <summary> /// <summary>
@ -295,10 +262,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -321,13 +285,27 @@ namespace MathNet.Numerics.Distributions
return Math.Log(Density(x)); return Math.Log(Density(x));
} }
#endregion
/// <summary>
/// Generates a sample from the Pareto distribution without doing parameter checking.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="scale">The scale parameter.</param>
/// <param name="shape">The shape parameter.</param>
/// <returns>a random number from the Pareto distribution.</returns>
internal static double SampleUnchecked(Random rnd, double scale, double shape)
{
return scale * Math.Pow(rnd.NextDouble(), -1.0 / shape);
}
/// <summary> /// <summary>
/// Draws a random sample from the distribution. /// Draws a random sample from the distribution.
/// </summary> /// </summary>
/// <returns>A random number from this distribution.</returns> /// <returns>A random number from this distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource); return SampleUnchecked(RandomSource, _scale, _shape);
} }
/// <summary> /// <summary>
@ -338,20 +316,45 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource); yield return SampleUnchecked(RandomSource, _scale, _shape);
} }
} }
#endregion /// <summary>
/// Generates a sample from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="scale">The scale parameter.</param>
/// <param name="shape">The shape parameter.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rnd, double scale, double shape)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(scale, shape))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, scale, shape);
}
/// <summary> /// <summary>
/// Generates a sample from the Pareto distribution without doing parameter checking. /// Generates a sequence of samples from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <returns>a random number from the Pareto distribution.</returns> /// <param name="scale">The scale parameter.</param>
private double DoSample(Random rnd) /// <param name="shape">The shape parameter.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rnd, double scale, double shape)
{ {
return _scale * Math.Pow(rnd.NextDouble(), -1.0 / _shape); if (Control.CheckDistributionParameters && !IsValidParameterSet(scale, shape))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleUnchecked(rnd, scale, shape);
}
} }
} }
} }

118
src/Numerics/Distributions/Continuous/Rayleigh.cs

@ -47,12 +47,12 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The scale parameter of the distribution. /// The scale parameter of the distribution.
/// </summary> /// </summary>
private double _scale; double _scale;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="Rayleigh"/> class. /// Initializes a new instance of the <see cref="Rayleigh"/> class.
@ -74,7 +74,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="scale">The scale parameter of the distribution.</param> /// <param name="scale">The scale parameter of the distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double scale) void SetParameters(double scale)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(scale)) if (Control.CheckDistributionParameters && !IsValidParameterSet(scale))
{ {
@ -89,7 +89,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
/// <param name="scale">The scale parameter of the distribution.</param> /// <param name="scale">The scale parameter of the distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double scale) static bool IsValidParameterSet(double scale)
{ {
if (scale <= 0) if (scale <= 0)
{ {
@ -109,15 +109,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return _scale; }
{
return _scale;
}
set set { SetParameters(value); }
{
SetParameters(value);
}
} }
/// <summary> /// <summary>
@ -136,10 +130,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -157,10 +148,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get get { return _scale * Math.Sqrt(Constants.PiOver2); }
{
return _scale * Math.Sqrt(Constants.PiOver2);
}
} }
/// <summary> /// <summary>
@ -168,10 +156,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get get { return (2.0 - Constants.PiOver2) * _scale * _scale; }
{
return (2.0 - Constants.PiOver2) * _scale * _scale;
}
} }
/// <summary> /// <summary>
@ -179,10 +164,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return Math.Sqrt(2.0 - Constants.PiOver2) * _scale; }
{
return Math.Sqrt(2.0 - Constants.PiOver2) * _scale;
}
} }
/// <summary> /// <summary>
@ -190,10 +172,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return 1.0 + Math.Log(_scale / Math.Sqrt(2)) + (Constants.EulerMascheroni / 2.0); }
{
return 1.0 + Math.Log(_scale / Math.Sqrt(2)) + (Constants.EulerMascheroni / 2.0);
}
} }
/// <summary> /// <summary>
@ -201,10 +180,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Skewness public double Skewness
{ {
get get { return (2.0 * Math.Sqrt(Constants.Pi) * (Constants.Pi - 3.0)) / Math.Pow(4.0 - Constants.Pi, 1.5); }
{
return (2.0 * Math.Sqrt(Constants.Pi) * (Constants.Pi - 3.0)) / Math.Pow(4.0 - Constants.Pi, 1.5);
}
} }
/// <summary> /// <summary>
@ -226,10 +202,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return _scale; }
{
return _scale;
}
} }
/// <summary> /// <summary>
@ -237,10 +210,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { return _scale * Math.Sqrt(Math.Log(4.0)); }
{
return _scale * Math.Sqrt(Math.Log(4.0));
}
} }
/// <summary> /// <summary>
@ -248,10 +218,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return 0.0; }
{
return 0.0;
}
} }
/// <summary> /// <summary>
@ -259,10 +226,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -285,13 +249,26 @@ namespace MathNet.Numerics.Distributions
return Math.Log(x / (_scale * _scale)) - (x * x / (2.0 * _scale * _scale)); return Math.Log(x / (_scale * _scale)) - (x * x / (2.0 * _scale * _scale));
} }
#endregion
/// <summary>
/// Generates a sample from the Rayleigh distribution without doing parameter checking.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="scale">The scale parameter.</param>
/// <returns>a random number from the Rayleigh distribution.</returns>
internal static double SampleUnchecked(Random rnd, double scale)
{
return scale * Math.Sqrt(-2.0 * Math.Log(rnd.NextDouble()));
}
/// <summary> /// <summary>
/// Draws a random sample from the distribution. /// Draws a random sample from the distribution.
/// </summary> /// </summary>
/// <returns>A random number from this distribution.</returns> /// <returns>A random number from this distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource); return SampleUnchecked(RandomSource, _scale);
} }
/// <summary> /// <summary>
@ -302,20 +279,43 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource); yield return SampleUnchecked(RandomSource, _scale);
} }
} }
#endregion /// <summary>
/// Generates a sample from the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="scale">The scale parameter.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rnd, double scale)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(scale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
return SampleUnchecked(rnd, scale);
}
/// <summary> /// <summary>
/// Generates a sample from the Rayleigh distribution without doing parameter checking. /// Generates a sequence of samples from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <returns>a random number from the Rayleigh distribution.</returns> /// <param name="scale">The scale parameter.</param>
private double DoSample(Random rnd) /// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rnd, double scale)
{ {
return _scale * Math.Sqrt(-2.0 * Math.Log(rnd.NextDouble())); if (Control.CheckDistributionParameters && !IsValidParameterSet(scale))
{
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
}
while (true)
{
yield return SampleUnchecked(rnd, scale);
}
} }
} }
} }

164
src/Numerics/Distributions/Continuous/Stable.cs

@ -47,27 +47,27 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// The stability parameter of the distribution. /// The stability parameter of the distribution.
/// </summary> /// </summary>
private double _alpha; double _alpha;
/// <summary> /// <summary>
/// The skewness parameter of the distribution. /// The skewness parameter of the distribution.
/// </summary> /// </summary>
private double _beta; double _beta;
/// <summary> /// <summary>
/// The scale parameter of the distribution. /// The scale parameter of the distribution.
/// </summary> /// </summary>
private double _scale; double _scale;
/// <summary> /// <summary>
/// The location parameter of the distribution. /// The location parameter of the distribution.
/// </summary> /// </summary>
private double _location; double _location;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the <see cref="Stable"/> class. /// Initializes a new instance of the <see cref="Stable"/> class.
@ -97,7 +97,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="beta">The skewness parameter of the distribution.</param> /// <param name="beta">The skewness parameter of the distribution.</param>
/// <param name="scale">The scale parameter of the distribution.</param> /// <param name="scale">The scale parameter of the distribution.</param>
/// <param name="location">The location parameter of the distribution.</param> /// <param name="location">The location parameter of the distribution.</param>
private void SetParameters(double alpha, double beta, double scale, double location) void SetParameters(double alpha, double beta, double scale, double location)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(alpha, beta, scale, location)) if (Control.CheckDistributionParameters && !IsValidParameterSet(alpha, beta, scale, location))
{ {
@ -118,7 +118,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="scale">The scale parameter of the distribution.</param> /// <param name="scale">The scale parameter of the distribution.</param>
/// <param name="location">The location parameter of the distribution.</param> /// <param name="location">The location parameter of the distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double alpha, double beta, double scale, double location) static bool IsValidParameterSet(double alpha, double beta, double scale, double location)
{ {
if (alpha <= 0 || alpha > 2) if (alpha <= 0 || alpha > 2)
{ {
@ -148,15 +148,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Alpha public double Alpha
{ {
get get { return _alpha; }
{
return _alpha;
}
set set { SetParameters(value, _beta, _scale, _location); }
{
SetParameters(value, _beta, _scale, _location);
}
} }
/// <summary> /// <summary>
@ -164,15 +158,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Beta public double Beta
{ {
get get { return _beta; }
{
return _beta;
}
set set { SetParameters(_alpha, value, _scale, _location); }
{
SetParameters(_alpha, value, _scale, _location);
}
} }
/// <summary> /// <summary>
@ -180,15 +168,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return _scale; }
{
return _scale;
}
set set { SetParameters(_alpha, _beta, value, _location); }
{
SetParameters(_alpha, _beta, value, _location);
}
} }
/// <summary> /// <summary>
@ -196,15 +178,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Location public double Location
{ {
get get { return _location; }
{
return _location;
}
set set { SetParameters(_alpha, _beta, _scale, value); }
{
SetParameters(_alpha, _beta, _scale, value);
}
} }
/// <summary> /// <summary>
@ -223,10 +199,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -293,10 +266,7 @@ namespace MathNet.Numerics.Distributions
/// <remarks>Always throws a not supported exception.</remarks> /// <remarks>Always throws a not supported exception.</remarks>
public double Entropy public double Entropy
{ {
get get { throw new NotSupportedException(); }
{
throw new NotSupportedException();
}
} }
/// <summary> /// <summary>
@ -351,7 +321,7 @@ namespace MathNet.Numerics.Distributions
/// <returns> /// <returns>
/// the cumulative density at <paramref name="x"/>. /// the cumulative density at <paramref name="x"/>.
/// </returns> /// </returns>
private static double LevyCumulativeDistribution(double scale, double location, double x) static double LevyCumulativeDistribution(double scale, double location, double x)
{ {
// The parameters scale and location must be correct // The parameters scale and location must be correct
return SpecialFunctions.Erfc(Math.Sqrt(scale / (2 * (x - location)))); return SpecialFunctions.Erfc(Math.Sqrt(scale / (2 * (x - location))));
@ -416,10 +386,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -454,7 +421,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="location">The location parameter of the distribution.</param> /// <param name="location">The location parameter of the distribution.</param>
/// <param name="x">The location at which to compute the density.</param> /// <param name="x">The location at which to compute the density.</param>
/// <returns>the density at <paramref name="x"/>.</returns> /// <returns>the density at <paramref name="x"/>.</returns>
private static double LevyDensity(double scale, double location, double x) static double LevyDensity(double scale, double location, double x)
{ {
// The parameters scale and location must be correct // The parameters scale and location must be correct
if (x < location) if (x < location)
@ -475,13 +442,51 @@ namespace MathNet.Numerics.Distributions
return Math.Log(Density(x)); return Math.Log(Density(x));
} }
#endregion
/// <summary>
/// Samples the distribution.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="alpha">The stability parameter of the distribution.</param>
/// <param name="beta">The skewness parameter of the distribution.</param>
/// <param name="scale">The scale parameter of the distribution.</param>
/// <param name="location">The location parameter of the distribution.</param>
/// <returns>a random number from the distribution.</returns>
internal static double SampleUnchecked(Random rnd, double alpha, double beta, double scale, double location)
{
var randTheta = ContinuousUniform.Sample(rnd, -Constants.PiOver2, Constants.PiOver2);
var randW = Exponential.Sample(rnd, 1.0);
if (!1.0.AlmostEqual(alpha))
{
var theta = (1.0 / alpha) * Math.Atan(beta * Math.Tan(Constants.PiOver2 * alpha));
var angle = alpha * (randTheta + theta);
var part1 = beta * Math.Tan(Constants.PiOver2 * alpha);
var factor = Math.Pow(1.0 + (part1 * part1), 1.0 / (2.0 * alpha));
var factor1 = Math.Sin(angle) / Math.Pow(Math.Cos(randTheta), (1.0 / alpha));
var factor2 = Math.Pow(Math.Cos(randTheta - angle) / randW, (1 - alpha) / alpha);
return location + scale * (factor * factor1 * factor2);
}
else
{
var part1 = Constants.PiOver2 + (beta * randTheta);
var summand = part1 * Math.Tan(randTheta);
var subtrahend = beta * Math.Log(Constants.PiOver2 * randW * Math.Cos(randTheta) / part1);
return location + scale * ((2.0 / Math.PI) * (summand - subtrahend));
}
}
/// <summary> /// <summary>
/// Draws a random sample from the distribution. /// Draws a random sample from the distribution.
/// </summary> /// </summary>
/// <returns>A random number from this distribution.</returns> /// <returns>A random number from this distribution.</returns>
public double Sample() public double Sample()
{ {
return DoSample(RandomSource, _alpha, _beta); return SampleUnchecked(RandomSource, _alpha, _beta, _scale, _location);
} }
/// <summary> /// <summary>
@ -492,43 +497,50 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return DoSample(RandomSource, _alpha, _beta); yield return SampleUnchecked(RandomSource, _alpha, _beta, _scale, _location);
} }
} }
/// <summary> /// <summary>
/// Samples the distribution. /// Generates a sample from the distribution.
/// </summary> /// </summary>
/// <param name="rnd">The random number generator to use.</param> /// <param name="rnd">The random number generator to use.</param>
/// <param name="alpha">The stability parameter of the distribution.</param> /// <param name="alpha">The stability parameter of the distribution.</param>
/// <param name="beta">The skewness parameter of the distribution.</param> /// <param name="beta">The skewness parameter of the distribution.</param>
/// <returns>a random number from the distribution.</returns> /// <param name="scale">The scale parameter of the distribution.</param>
private static double DoSample(Random rnd, double alpha, double beta) /// <param name="location">The location parameter of the distribution.</param>
/// <returns>a sample from the distribution.</returns>
public static double Sample(Random rnd, double alpha, double beta, double scale, double location)
{ {
var randTheta = ContinuousUniform.Sample(rnd, -Constants.PiOver2, Constants.PiOver2); if (Control.CheckDistributionParameters && !IsValidParameterSet(alpha, beta, scale, location))
var randW = Exponential.Sample(rnd, 1.0);
if (!1.0.AlmostEqual(alpha))
{ {
var theta = (1.0 / alpha) * Math.Atan(beta * Math.Tan(Constants.PiOver2 * alpha)); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
var angle = alpha * (randTheta + theta); }
var part1 = beta * Math.Tan(Constants.PiOver2 * alpha);
var factor = Math.Pow(1.0 + (part1 * part1), 1.0 / (2.0 * alpha)); return SampleUnchecked(rnd, location, scale, scale, location);
var factor1 = Math.Sin(angle) / Math.Pow(Math.Cos(randTheta), (1.0 / alpha)); }
var factor2 = Math.Pow(Math.Cos(randTheta - angle) / randW, (1 - alpha) / alpha);
return factor * factor1 * factor2; /// <summary>
} /// Generates a sequence of samples from the distribution.
else /// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="alpha">The stability parameter of the distribution.</param>
/// <param name="beta">The skewness parameter of the distribution.</param>
/// <param name="scale">The scale parameter of the distribution.</param>
/// <param name="location">The location parameter of the distribution.</param>
/// <returns>a sequence of samples from the distribution.</returns>
public static IEnumerable<double> Samples(Random rnd, double alpha, double beta, double scale, double location)
{
if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale, scale, location))
{ {
var part1 = Constants.PiOver2 + (beta * randTheta); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
var summand = part1 * Math.Tan(randTheta); }
var subtrahend = beta * Math.Log(Constants.PiOver2 * randW * Math.Cos(randTheta) / part1);
return (2.0 / Math.PI) * (summand - subtrahend); while (true)
{
yield return SampleUnchecked(rnd, location, scale, scale, location);
} }
} }
#endregion
} }
} }

131
src/Numerics/Distributions/Continuous/StudentT.cs

@ -55,29 +55,30 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Keeps track of the location of the Student t-distribution. /// Keeps track of the location of the Student t-distribution.
/// </summary> /// </summary>
private double _location; double _location;
/// <summary> /// <summary>
/// Keeps track of the degrees of freedom for the Student t-distribution. /// Keeps track of the degrees of freedom for the Student t-distribution.
/// </summary> /// </summary>
private double _dof; double _dof;
/// <summary> /// <summary>
/// Keeps track of the scale for the Student t-distribution. /// Keeps track of the scale for the Student t-distribution.
/// </summary> /// </summary>
private double _scale; double _scale;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the StudentT class. This is a Student t-distribution with location 0.0 /// Initializes a new instance of the StudentT class. This is a Student t-distribution with location 0.0
/// scale 1.0 and degrees of freedom 1. The distribution will /// scale 1.0 and degrees of freedom 1. The distribution will
/// be initialized with the default <seealso cref="System.Random"/> random number generator. /// be initialized with the default <seealso cref="System.Random"/> random number generator.
/// </summary> /// </summary>
public StudentT() : this(0.0, 1.0, 1.0) public StudentT()
: this(0.0, 1.0, 1.0)
{ {
} }
@ -111,7 +112,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="scale">The scale of the Student t-distribution.</param> /// <param name="scale">The scale of the Student t-distribution.</param>
/// <param name="dof">The degrees of freedom for the Student t-distribution.</param> /// <param name="dof">The degrees of freedom for the Student t-distribution.</param>
/// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters are valid, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double location, double scale, double dof) static bool IsValidParameterSet(double location, double scale, double dof)
{ {
if (scale <= 0.0 || dof <= 0.0 || Double.IsNaN(scale) || Double.IsNaN(location) || Double.IsNaN(dof)) if (scale <= 0.0 || dof <= 0.0 || Double.IsNaN(scale) || Double.IsNaN(location) || Double.IsNaN(dof))
{ {
@ -128,7 +129,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="scale">The scale of the Student t-distribution.</param> /// <param name="scale">The scale of the Student t-distribution.</param>
/// <param name="dof">The degrees of freedom for the Student t-distribution.</param> /// <param name="dof">The degrees of freedom for the Student t-distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double location, double scale, double dof) void SetParameters(double location, double scale, double dof)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale, dof)) if (Control.CheckDistributionParameters && !IsValidParameterSet(location, scale, dof))
{ {
@ -145,15 +146,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Location public double Location
{ {
get get { return _location; }
{
return _location;
}
set set { SetParameters(value, _scale, _dof); }
{
SetParameters(value, _scale, _dof);
}
} }
/// <summary> /// <summary>
@ -161,15 +156,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return _scale; }
{
return _scale;
}
set set { SetParameters(_location, value, _dof); }
{
SetParameters(_location, value, _dof);
}
} }
/// <summary> /// <summary>
@ -177,15 +166,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double DegreesOfFreedom public double DegreesOfFreedom
{ {
get get { return _dof; }
{
return _dof;
}
set set { SetParameters(_location, _scale, value); }
{
SetParameters(_location, _scale, value);
}
} }
#region IDistribution implementation #region IDistribution implementation
@ -195,10 +178,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -216,10 +196,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get get { return _dof > 1.0 ? _location : Double.NaN; }
{
return _dof > 1.0 ? _location : Double.NaN;
}
} }
/// <summary> /// <summary>
@ -305,10 +282,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mode public double Mode
{ {
get get { return _location; }
{
return _location;
}
} }
/// <summary> /// <summary>
@ -316,10 +290,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { return _location; }
{
return _location;
}
} }
/// <summary> /// <summary>
@ -327,10 +298,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Minimum public double Minimum
{ {
get get { return Double.NegativeInfinity; }
{
return Double.NegativeInfinity;
}
} }
/// <summary> /// <summary>
@ -338,10 +306,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Maximum public double Maximum
{ {
get get { return Double.PositiveInfinity; }
{
return Double.PositiveInfinity;
}
} }
/// <summary> /// <summary>
@ -359,9 +324,9 @@ namespace MathNet.Numerics.Distributions
var d = (x - _location) / _scale; var d = (x - _location) / _scale;
return Math.Exp(SpecialFunctions.GammaLn((_dof + 1.0) / 2.0) - SpecialFunctions.GammaLn(_dof / 2.0)) return Math.Exp(SpecialFunctions.GammaLn((_dof + 1.0) / 2.0) - SpecialFunctions.GammaLn(_dof / 2.0))
* Math.Pow(1.0 + (d * d / _dof), -0.5 * (_dof + 1.0)) * Math.Pow(1.0 + (d * d / _dof), -0.5 * (_dof + 1.0))
/ Math.Sqrt(_dof * Math.PI) / Math.Sqrt(_dof * Math.PI)
/ _scale; / _scale;
} }
/// <summary> /// <summary>
@ -379,9 +344,9 @@ namespace MathNet.Numerics.Distributions
var d = (x - _location) / _scale; var d = (x - _location) / _scale;
return SpecialFunctions.GammaLn((_dof + 1.0) / 2.0) return SpecialFunctions.GammaLn((_dof + 1.0) / 2.0)
- (0.5 * ((_dof + 1.0) * Math.Log(1.0 + (d * d / _dof)))) - (0.5 * ((_dof + 1.0) * Math.Log(1.0 + (d * d / _dof))))
- SpecialFunctions.GammaLn(_dof / 2.0) - SpecialFunctions.GammaLn(_dof / 2.0)
- (0.5 * Math.Log(_dof * Math.PI)) - Math.Log(_scale); - (0.5 * Math.Log(_dof * Math.PI)) - Math.Log(_scale);
} }
/// <summary> /// <summary>
@ -403,13 +368,32 @@ namespace MathNet.Numerics.Distributions
return x <= _location ? ib : 1.0 - ib; return x <= _location ? ib : 1.0 - ib;
} }
#endregion
/// <summary>
/// Samples student-t distributed random variables.
/// </summary>
/// <remarks>The algorithm is method 2 in section 5, chapter 9
/// in L. Devroye's "Non-Uniform Random Variate Generation"</remarks>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="location">The location of the Student t-distribution.</param>
/// <param name="scale">The scale of the Student t-distribution.</param>
/// <param name="dof">The degrees of freedom for the standard student-t distribution.</param>
/// <returns>a random number from the standard student-t distribution.</returns>
internal static double SampleUnchecked(Random rnd, double location, double scale, double dof)
{
var n = Normal.SampleUncheckedBoxMuller(rnd).Item1;
var g = Gamma.SampleUnchecked(rnd, 0.5 * dof, 0.5);
return location + (scale * n * Math.Sqrt(dof / g));
}
/// <summary> /// <summary>
/// Generates a sample from the Student t-distribution. /// Generates a sample from the Student t-distribution.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return _location + (_scale * Sample(RandomSource, _dof)); return SampleUnchecked(RandomSource, _location, _scale, _dof);
} }
/// <summary> /// <summary>
@ -420,12 +404,10 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return _location + (_scale * Sample(RandomSource, _dof)); yield return SampleUnchecked(RandomSource, _location, _scale, _dof);
} }
} }
#endregion
/// <summary> /// <summary>
/// Generates a sample from the Student t-distribution. /// Generates a sample from the Student t-distribution.
/// </summary> /// </summary>
@ -441,7 +423,7 @@ namespace MathNet.Numerics.Distributions
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return location + (scale * Sample(rng, dof)); return SampleUnchecked(rng, location, scale, dof);
} }
/// <summary> /// <summary>
@ -461,23 +443,8 @@ namespace MathNet.Numerics.Distributions
while (true) while (true)
{ {
yield return location + (scale * Sample(rng, dof)); yield return SampleUnchecked(rng, location, scale, dof);
} }
} }
/// <summary>
/// Samples standard student-t distributed random variables.
/// </summary>
/// <remarks>The algorithm is method 2 in section 5, chapter 9
/// in L. Devroye's "Non-Uniform Random Variate Generation"</remarks>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="dof">The degrees of freedom for the standard student-t distribution.</param>
/// <returns>a random number from the standard student-t distribution.</returns>
internal static double Sample(Random rnd, double dof)
{
var n = Normal.SampleBoxMuller(rnd).Item1;
var g = Gamma.Sample(rnd, 0.5 * dof, 0.5);
return Math.Sqrt(dof / g) * n;
}
} }
} }

103
src/Numerics/Distributions/Continuous/Weibull.cs

@ -50,12 +50,12 @@ namespace MathNet.Numerics.Distributions
/// <summary> /// <summary>
/// Weibull shape parameter. /// Weibull shape parameter.
/// </summary> /// </summary>
private double _shape; double _shape;
/// <summary> /// <summary>
/// Weibull inverse scale parameter. /// Weibull inverse scale parameter.
/// </summary> /// </summary>
private double _scale; double _scale;
/// <summary> /// <summary>
/// Reusable intermediate result 1 / (<see cref="_scale"/> ^ <see cref="_shape"/>) /// Reusable intermediate result 1 / (<see cref="_scale"/> ^ <see cref="_shape"/>)
@ -64,12 +64,12 @@ namespace MathNet.Numerics.Distributions
/// By caching this parameter we can get slightly better numerics precision /// By caching this parameter we can get slightly better numerics precision
/// in certain constellations without any additional computations. /// in certain constellations without any additional computations.
/// </remarks> /// </remarks>
private double _scalePowShapeInv; double _scalePowShapeInv;
/// <summary> /// <summary>
/// The distribution's random number generator. /// The distribution's random number generator.
/// </summary> /// </summary>
private Random _random; Random _random;
/// <summary> /// <summary>
/// Initializes a new instance of the Weibull class. /// Initializes a new instance of the Weibull class.
@ -97,7 +97,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="shape">The shape of the Weibull distribution.</param> /// <param name="shape">The shape of the Weibull distribution.</param>
/// <param name="scale">The scale of the Weibull distribution.</param> /// <param name="scale">The scale of the Weibull distribution.</param>
/// <returns><c>true</c> when the parameters positive valid floating point numbers, <c>false</c> otherwise.</returns> /// <returns><c>true</c> when the parameters positive valid floating point numbers, <c>false</c> otherwise.</returns>
private static bool IsValidParameterSet(double shape, double scale) static bool IsValidParameterSet(double shape, double scale)
{ {
if (shape <= 0.0 || scale <= 0.0 || Double.IsNaN(shape) || Double.IsNaN(scale)) if (shape <= 0.0 || scale <= 0.0 || Double.IsNaN(shape) || Double.IsNaN(scale))
{ {
@ -113,7 +113,7 @@ namespace MathNet.Numerics.Distributions
/// <param name="shape">The shape of the Weibull distribution.</param> /// <param name="shape">The shape of the Weibull distribution.</param>
/// <param name="scale">The inverse scale of the Weibull distribution.</param> /// <param name="scale">The inverse scale of the Weibull distribution.</param>
/// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception> /// <exception cref="ArgumentOutOfRangeException">When the parameters don't pass the <see cref="IsValidParameterSet"/> function.</exception>
private void SetParameters(double shape, double scale) void SetParameters(double shape, double scale)
{ {
if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, scale)) if (Control.CheckDistributionParameters && !IsValidParameterSet(shape, scale))
{ {
@ -130,15 +130,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Shape public double Shape
{ {
get get { return _shape; }
{
return _shape;
}
set set { SetParameters(value, _scale); }
{
SetParameters(value, _scale);
}
} }
/// <summary> /// <summary>
@ -146,15 +140,9 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Scale public double Scale
{ {
get get { return _scale; }
{
return _scale;
}
set set { SetParameters(_shape, value); }
{
SetParameters(_shape, value);
}
} }
#region IDistribution implementation #region IDistribution implementation
@ -164,10 +152,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public Random RandomSource public Random RandomSource
{ {
get get { return _random; }
{
return _random;
}
set set
{ {
@ -185,10 +170,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Mean public double Mean
{ {
get get { return _scale * SpecialFunctions.Gamma(1.0 + (1.0 / _shape)); }
{
return _scale * SpecialFunctions.Gamma(1.0 + (1.0 / _shape));
}
} }
/// <summary> /// <summary>
@ -196,10 +178,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Variance public double Variance
{ {
get get { return (_scale * _scale * SpecialFunctions.Gamma(1.0 + (2.0 / _shape))) - (Mean * Mean); }
{
return (_scale * _scale * SpecialFunctions.Gamma(1.0 + (2.0 / _shape))) - (Mean * Mean);
}
} }
/// <summary> /// <summary>
@ -207,10 +186,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double StdDev public double StdDev
{ {
get get { return Math.Sqrt(Variance); }
{
return Math.Sqrt(Variance);
}
} }
/// <summary> /// <summary>
@ -218,10 +194,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Entropy public double Entropy
{ {
get get { return (Constants.EulerMascheroni * (1.0 - (1.0 / _shape))) + Math.Log(_scale / _shape) + 1.0; }
{
return (Constants.EulerMascheroni * (1.0 - (1.0 / _shape))) + Math.Log(_scale / _shape) + 1.0;
}
} }
/// <summary> /// <summary>
@ -238,6 +211,7 @@ namespace MathNet.Numerics.Distributions
return ((_scale * _scale * _scale * SpecialFunctions.Gamma(1.0 + (3.0 / _shape))) - (3.0 * sigma2 * mu) - (mu * mu * mu)) / sigma3; return ((_scale * _scale * _scale * SpecialFunctions.Gamma(1.0 + (3.0 / _shape))) - (3.0 * sigma2 * mu) - (mu * mu * mu)) / sigma3;
} }
} }
#endregion #endregion
#region IContinuousDistribution implementation #region IContinuousDistribution implementation
@ -263,10 +237,7 @@ namespace MathNet.Numerics.Distributions
/// </summary> /// </summary>
public double Median public double Median
{ {
get get { return _scale * Math.Pow(Constants.Ln2, 1.0 / _shape); }
{
return _scale * Math.Pow(Constants.Ln2, 1.0 / _shape);
}
} }
/// <summary> /// <summary>
@ -340,13 +311,29 @@ namespace MathNet.Numerics.Distributions
return -SpecialFunctions.ExponentialMinusOne(-Math.Pow(x, _shape) * _scalePowShapeInv); return -SpecialFunctions.ExponentialMinusOne(-Math.Pow(x, _shape) * _scalePowShapeInv);
} }
#endregion
/// <summary>
/// Generates one sample from the Weibull distribution. This method doesn't perform
/// any parameter checks.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="shape">The shape of the Weibull distribution.</param>
/// <param name="scale">The scale of the Weibull distribution.</param>
/// <returns>A sample from a Weibull distributed random variable.</returns>
internal static double SampleUnchecked(Random rnd, double shape, double scale)
{
var x = rnd.NextDouble();
return scale * Math.Pow(-Math.Log(x), 1.0 / shape);
}
/// <summary> /// <summary>
/// Generates a sample from the Weibull distribution. /// Generates a sample from the Weibull distribution.
/// </summary> /// </summary>
/// <returns>a sample from the distribution.</returns> /// <returns>a sample from the distribution.</returns>
public double Sample() public double Sample()
{ {
return SampleWeibull(RandomSource, _shape, _scale); return SampleUnchecked(RandomSource, _shape, _scale);
} }
/// <summary> /// <summary>
@ -357,12 +344,10 @@ namespace MathNet.Numerics.Distributions
{ {
while (true) while (true)
{ {
yield return SampleWeibull(RandomSource, _shape, _scale); yield return SampleUnchecked(RandomSource, _shape, _scale);
} }
} }
#endregion
/// <summary> /// <summary>
/// Generates a sample from the Weibull distribution. /// Generates a sample from the Weibull distribution.
/// </summary> /// </summary>
@ -377,7 +362,7 @@ namespace MathNet.Numerics.Distributions
throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
} }
return SampleWeibull(rng, shape, scale); return SampleUnchecked(rng, shape, scale);
} }
/// <summary> /// <summary>
@ -396,22 +381,8 @@ namespace MathNet.Numerics.Distributions
while (true) while (true)
{ {
yield return SampleWeibull(rng, shape, scale); yield return SampleUnchecked(rng, shape, scale);
} }
} }
/// <summary>
/// Generates one sample from the Weibull distribution. This method doesn't perform
/// any parameter checks.
/// </summary>
/// <param name="rnd">The random number generator to use.</param>
/// <param name="shape">The shape of the Weibull distribution.</param>
/// <param name="scale">The scale of the Weibull distribution.</param>
/// <returns>A sample from a Weibull distributed random variable.</returns>
internal static double SampleWeibull(Random rnd, double shape, double scale)
{
var x = rnd.NextDouble();
return scale * Math.Pow(-Math.Log(x), 1.0 / shape);
}
} }
} }

2
src/Numerics/Distributions/Multivariate/MatrixNormal.cs

@ -326,7 +326,7 @@ namespace MathNet.Numerics.Distributions
var v = new DenseVector(count, 0.0); var v = new DenseVector(count, 0.0);
for (var d = 0; d < count; d += 2) for (var d = 0; d < count; d += 2)
{ {
var sample = Normal.SampleBoxMuller(rnd); var sample = Normal.SampleUncheckedBoxMuller(rnd);
v[d] = sample.Item1; v[d] = sample.Item1;
if (d + 1 < count) if (d + 1 < count)
{ {

Loading…
Cancel
Save