diff --git a/src/Numerics/Distributions/Cauchy.cs b/src/Numerics/Distributions/Cauchy.cs index bea65ade..3435834f 100644 --- a/src/Numerics/Distributions/Cauchy.cs +++ b/src/Numerics/Distributions/Cauchy.cs @@ -258,7 +258,7 @@ namespace MathNet.Numerics.Distributions /// A random number from this distribution. public double Sample() { - return SampleUnchecked(_random, _location, _scale); + return _location + _scale*Math.Tan(Constants.Pi*(_random.NextDouble() - 0.5)); } /// @@ -269,23 +269,10 @@ namespace MathNet.Numerics.Distributions { while (true) { - yield return SampleUnchecked(_random, _location, _scale); + yield return _location + _scale*Math.Tan(Constants.Pi*(_random.NextDouble() - 0.5)); } } - /// - /// Samples the distribution. - /// - /// The random number generator to use. - /// The location (x0) of the distribution. - /// The scale (γ) of the distribution. Range: γ > 0. - /// a random number from the distribution. - static double SampleUnchecked(System.Random rnd, double location, double scale) - { - var u = rnd.NextDouble(); - return location + (scale*Math.Tan(Constants.Pi*(u - 0.5))); - } - /// /// Computes the probability density of the distribution (PDF) at x, i.e. ∂P(X ≤ x)/∂x. /// @@ -359,7 +346,7 @@ namespace MathNet.Numerics.Distributions { if (scale <= 0.0) throw new ArgumentOutOfRangeException("scale", Resources.InvalidDistributionParameters); - return SampleUnchecked(rnd, location, scale); + return location + scale*Math.Tan(Constants.Pi*(rnd.NextDouble() - 0.5)); } /// @@ -375,7 +362,7 @@ namespace MathNet.Numerics.Distributions while (true) { - yield return SampleUnchecked(rnd, location, scale); + yield return location + scale*Math.Tan(Constants.Pi*(rnd.NextDouble() - 0.5)); } } } diff --git a/src/Numerics/Distributions/ContinuousUniform.cs b/src/Numerics/Distributions/ContinuousUniform.cs index 0a729854..2cab6128 100644 --- a/src/Numerics/Distributions/ContinuousUniform.cs +++ b/src/Numerics/Distributions/ContinuousUniform.cs @@ -92,17 +92,6 @@ namespace MathNet.Numerics.Distributions return "ContinuousUniform(Lower = " + _lower + ", Upper = " + _upper + ")"; } - /// - /// Checks whether the parameters of the distribution are valid. - /// - /// Lower bound. Range: lower ≤ upper. - /// Upper bound. Range: lower ≤ upper. - /// true when the parameters are valid, false otherwise. - static bool IsValidParameterSet(double lower, double upper) - { - return upper >= lower && !Double.IsNaN(upper) && !Double.IsNaN(lower); - } - /// /// Sets the parameters of the distribution after checking their validity. /// @@ -111,7 +100,7 @@ namespace MathNet.Numerics.Distributions /// When the parameters are out of range. void SetParameters(double lower, double upper) { - if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper)) + if (upper < lower || Double.IsNaN(upper) || Double.IsNaN(lower)) { throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); } @@ -227,14 +216,10 @@ namespace MathNet.Numerics.Distributions /// /// The location at which to compute the density. /// the density at . + /// public double Density(double x) { - if (x >= _lower && x <= _upper) - { - return 1.0/(_upper - _lower); - } - - return 0.0; + return x < _lower || x > _upper ? 0.0 : 1.0/(_upper - _lower); } /// @@ -242,14 +227,10 @@ namespace MathNet.Numerics.Distributions /// /// The location at which to compute the log density. /// the log density at . + /// public double DensityLn(double x) { - if (x >= _lower && x <= _upper) - { - return -Math.Log(_upper - _lower); - } - - return Double.NegativeInfinity; + return x < _lower || x > _upper ? Double.NegativeInfinity : -Math.Log(_upper - _lower); } /// @@ -257,31 +238,22 @@ namespace MathNet.Numerics.Distributions /// /// The location at which to compute the cumulative distribution function. /// the cumulative distribution at location . + /// public double CumulativeDistribution(double x) { - if (x <= _lower) - { - return 0.0; - } - - if (x >= _upper) - { - return 1.0; - } - - return (x - _lower)/(_upper - _lower); + return x <= _lower ? 0.0 : x >= _upper ? 1.0 : (x - _lower)/(_upper - _lower); } /// - /// Generates one sample from the ContinuousUniform distribution without parameter checking. + /// Computes the inverse of the cumulative distribution function (InvCDF) for the distribution + /// at the given probability. This is also known as the quantile or percent point function. /// - /// The random number generator to use. - /// Lower bound. Range: lower ≤ upper. - /// Upper bound. Range: lower ≤ upper. - /// a uniformly distributed random number. - static double SampleUnchecked(System.Random rnd, double lower, double upper) + /// The location at which to compute the inverse cumulative density. + /// the inverse cumulative density at . + /// + public double InverseCumulativeDistribution(double p) { - return lower + (rnd.NextDouble()*(upper - lower)); + return p <= 0.0 ? _lower : p >= 1.0 ? _upper : _lower*(1.0 - p) + _upper*p; } /// @@ -290,7 +262,7 @@ namespace MathNet.Numerics.Distributions /// a sample from the distribution. public double Sample() { - return SampleUnchecked(_random, _lower, _upper); + return _lower + _random.NextDouble()*(_upper - _lower); } /// @@ -301,10 +273,71 @@ namespace MathNet.Numerics.Distributions { while (true) { - yield return SampleUnchecked(_random, _lower, _upper); + yield return _lower + _random.NextDouble()*(_upper - _lower); } } + /// + /// Computes the probability density of the distribution (PDF) at x, i.e. ∂P(X ≤ x)/∂x. + /// + /// Lower bound. Range: lower ≤ upper. + /// Upper bound. Range: lower ≤ upper. + /// The location at which to compute the density. + /// the density at . + /// + public static double PDF(double lower, double upper, double x) + { + if (upper < lower) throw new ArgumentOutOfRangeException("upper", Resources.InvalidDistributionParameters); + + return x < lower || x > upper ? 0.0 : 1.0/(upper - lower); + } + + /// + /// Computes the log probability density of the distribution (lnPDF) at x, i.e. ln(∂P(X ≤ x)/∂x). + /// + /// Lower bound. Range: lower ≤ upper. + /// Upper bound. Range: lower ≤ upper. + /// The location at which to compute the density. + /// the log density at . + /// + public static double PDFLn(double lower, double upper, double x) + { + if (upper < lower) throw new ArgumentOutOfRangeException("upper", Resources.InvalidDistributionParameters); + + return x < lower || x > upper ? Double.NegativeInfinity : -Math.Log(upper - lower); + } + + /// + /// Computes the cumulative distribution (CDF) of the distribution at x, i.e. P(X ≤ x). + /// + /// The location at which to compute the cumulative distribution function. + /// Lower bound. Range: lower ≤ upper. + /// Upper bound. Range: lower ≤ upper. + /// the cumulative distribution at location . + /// + public static double CDF(double lower, double upper, double x) + { + if (upper < lower) throw new ArgumentOutOfRangeException("upper", Resources.InvalidDistributionParameters); + + return x <= lower ? 0.0 : x >= upper ? 1.0 : (x - lower)/(upper - lower); + } + + /// + /// Computes the inverse of the cumulative distribution function (InvCDF) for the distribution + /// at the given probability. This is also known as the quantile or percent point function. + /// + /// The location at which to compute the inverse cumulative density. + /// Lower bound. Range: lower ≤ upper. + /// Upper bound. Range: lower ≤ upper. + /// the inverse cumulative density at . + /// + public static double InvCDF(double lower, double upper, double p) + { + if (upper < lower) throw new ArgumentOutOfRangeException("upper", Resources.InvalidDistributionParameters); + + return p <= 0.0 ? lower : p >= 1.0 ? upper : lower*(1.0 - p) + upper*p; + } + /// /// Generates a sample from the ContinuousUniform distribution. /// @@ -314,12 +347,9 @@ namespace MathNet.Numerics.Distributions /// a uniformly distributed sample. public static double Sample(System.Random rnd, double lower, double upper) { - if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper)) - { - throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); - } + if (upper < lower) throw new ArgumentOutOfRangeException("upper", Resources.InvalidDistributionParameters); - return SampleUnchecked(rnd, lower, upper); + return lower + rnd.NextDouble()*(upper - lower); } /// @@ -331,14 +361,11 @@ namespace MathNet.Numerics.Distributions /// a sequence of uniformly distributed samples. public static IEnumerable Samples(System.Random rnd, double lower, double upper) { - if (Control.CheckDistributionParameters && !IsValidParameterSet(lower, upper)) - { - throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters); - } + if (upper < lower) throw new ArgumentOutOfRangeException("upper", Resources.InvalidDistributionParameters); while (true) { - yield return SampleUnchecked(rnd, lower, upper); + yield return lower + rnd.NextDouble()*(upper - lower); } } } diff --git a/src/UnitTests/DistributionTests/Continuous/ContinuousUniformTests.cs b/src/UnitTests/DistributionTests/Continuous/ContinuousUniformTests.cs index 363a47d5..1b1bac88 100644 --- a/src/UnitTests/DistributionTests/Continuous/ContinuousUniformTests.cs +++ b/src/UnitTests/DistributionTests/Continuous/ContinuousUniformTests.cs @@ -275,11 +275,13 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous var x = i - 5.0; if (x >= lower && x <= upper) { - Assert.AreEqual(1.0 / (upper - lower), n.Density(x)); + Assert.AreEqual(1.0/(upper - lower), n.Density(x)); + Assert.AreEqual(1.0/(upper - lower), ContinuousUniform.PDF(lower, upper, x)); } else { Assert.AreEqual(0.0, n.Density(x)); + Assert.AreEqual(0.0, ContinuousUniform.PDF(lower, upper, x)); } } } @@ -304,10 +306,12 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous if (x >= lower && x <= upper) { Assert.AreEqual(-Math.Log(upper - lower), n.DensityLn(x)); + Assert.AreEqual(-Math.Log(upper - lower), ContinuousUniform.PDFLn(lower, upper, x)); } else { Assert.AreEqual(double.NegativeInfinity, n.DensityLn(x)); + Assert.AreEqual(double.NegativeInfinity, ContinuousUniform.PDFLn(lower, upper, x)); } } } @@ -390,14 +394,51 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous if (x <= lower) { Assert.AreEqual(0.0, n.CumulativeDistribution(x)); + Assert.AreEqual(0.0, ContinuousUniform.CDF(lower, upper, x)); } else if (x >= upper) { Assert.AreEqual(1.0, n.CumulativeDistribution(x)); + Assert.AreEqual(1.0, ContinuousUniform.CDF(lower, upper, x)); } else { - Assert.AreEqual((x - lower) / (upper - lower), n.CumulativeDistribution(x)); + Assert.AreEqual((x - lower)/(upper - lower), n.CumulativeDistribution(x)); + Assert.AreEqual((x - lower)/(upper - lower), ContinuousUniform.CDF(lower, upper, x)); + } + } + } + + /// + /// Validate inverse cumulative distribution. + /// + /// Lower bound. + /// Upper bound. + [TestCase(0.0, 0.0)] + [TestCase(0.0, 0.1)] + [TestCase(0.0, 1.0)] + [TestCase(0.0, 10.0)] + [TestCase(-5.0, 100.0)] + public void ValidateInverseCumulativeDistribution(double lower, double upper) + { + var n = new ContinuousUniform(lower, upper); + for (var i = 0; i < 11; i++) + { + var x = i - 5.0; + if (x <= lower) + { + Assert.AreEqual(lower, n.InverseCumulativeDistribution(0.0), 1e-12); + Assert.AreEqual(lower, ContinuousUniform.InvCDF(lower, upper, 0.0), 1e-12); + } + else if (x >= upper) + { + Assert.AreEqual(upper, n.InverseCumulativeDistribution(1.0), 1e-12); + Assert.AreEqual(upper, ContinuousUniform.InvCDF(lower, upper, 1.0), 1e-12); + } + else + { + Assert.AreEqual(x, n.InverseCumulativeDistribution((x - lower)/(upper - lower)), 1e-12); + Assert.AreEqual(x, ContinuousUniform.InvCDF(lower, upper, (x - lower)/(upper - lower)), 1e-12); } } }