diff --git a/src/Numerics/Distributions/StudentT.cs b/src/Numerics/Distributions/StudentT.cs
index c52f4e36..4d4b50ce 100644
--- a/src/Numerics/Distributions/StudentT.cs
+++ b/src/Numerics/Distributions/StudentT.cs
@@ -32,6 +32,7 @@ using System;
using System.Collections.Generic;
using MathNet.Numerics.Properties;
using MathNet.Numerics.Random;
+using MathNet.Numerics.RootFinding;
namespace MathNet.Numerics.Distributions
{
@@ -304,6 +305,19 @@ namespace MathNet.Numerics.Distributions
return CDF(_location, _scale, _freedom, x);
}
+ ///
+ /// 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.
+ /// the inverse cumulative density at .
+ ///
+ /// WARNING: currently not an explicit implementation, hence slow and unreliable.
+ public double InverseCumulativeDistribution(double p)
+ {
+ return InvCDF(_location, _scale, _freedom, p);
+ }
+
///
/// Samples student-t distributed random variables.
///
@@ -409,6 +423,35 @@ namespace MathNet.Numerics.Distributions
return x <= location ? ib : 1.0 - ib;
}
+ ///
+ /// 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.
+ /// The location (μ) of the distribution.
+ /// The scale (σ) of the distribution. Range: σ > 0.
+ /// The degrees of freedom (ν) for the distribution. Range: ν > 0.
+ /// the inverse cumulative density at .
+ ///
+ /// WARNING: currently not an explicit implementation, hence slow and unreliable.
+ public static double InvCDF(double location, double scale, double freedom, double p)
+ {
+ if (scale <= 0.0 || freedom <= 0.0) throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
+
+ // TODO JVG we can probably do a better job for Cauchy special case
+ if (Double.IsPositiveInfinity(freedom)) return Normal.InvCDF(location, scale, p);
+ if (p == 0.5d) return location;
+
+ // TODO PERF: We must implement this explicitly instead of solving for CDF^-1
+ return Brent.FindRoot(x =>
+ {
+ var k = (x - location)/scale;
+ var h = freedom/(freedom + (k*k));
+ var ib = 0.5*SpecialFunctions.BetaRegularized(freedom/2.0, 0.5, h);
+ return x <= location ? ib - p : 1.0 - ib - p;
+ }, -800, 800, accuracy: 1e-8);
+ }
+
///
/// Generates a sample from the Student t-distribution.
///
diff --git a/src/UnitTests/DistributionTests/Continuous/NormalTests.cs b/src/UnitTests/DistributionTests/Continuous/NormalTests.cs
index 30e4d57f..bf163821 100644
--- a/src/UnitTests/DistributionTests/Continuous/NormalTests.cs
+++ b/src/UnitTests/DistributionTests/Continuous/NormalTests.cs
@@ -458,7 +458,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
/// Validate cumulative distribution.
///
/// Input X value.
- /// Expected value.
+ /// Expected value.
[TestCase(Double.NegativeInfinity, 0.0)]
[TestCase(-5.0, 0.00000028665157187919391167375233287464535385442301361187883)]
[TestCase(-2.0, 0.0002326290790355250363499258867279847735487493358890356)]
@@ -469,18 +469,18 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
[TestCase(6.0, 0.69146246127401310363770461060833773988360217555457859)]
[TestCase(10.0, 0.9937903346742238648330218954258077788721022530769078)]
[TestCase(Double.PositiveInfinity, 1.0)]
- public void ValidateCumulativeDistribution(double x, double f)
+ public void ValidateCumulativeDistribution(double x, double p)
{
var n = Normal.WithMeanStdDev(5.0, 2.0);
- AssertHelpers.AlmostEqualRelative(f, n.CumulativeDistribution(x), 9);
- AssertHelpers.AlmostEqualRelative(f, Normal.CDF(5.0, 2.0, x), 9);
+ AssertHelpers.AlmostEqualRelative(p, n.CumulativeDistribution(x), 9);
+ AssertHelpers.AlmostEqualRelative(p, Normal.CDF(5.0, 2.0, x), 9);
}
///
/// Validate inverse cumulative distribution.
///
/// Input X value.
- /// Expected value.
+ /// Expected value.
[TestCase(Double.NegativeInfinity, 0.0)]
[TestCase(-5.0, 0.00000028665157187919391167375233287464535385442301361187883)]
[TestCase(-2.0, 0.0002326290790355250363499258867279847735487493358890356)]
@@ -491,11 +491,11 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
[TestCase(6.0, .69146246127401310363770461060833773988360217555457859)]
[TestCase(10.0, 0.9937903346742238648330218954258077788721022530769078)]
[TestCase(Double.PositiveInfinity, 1.0)]
- public void ValidateInverseCumulativeDistribution(double x, double f)
+ public void ValidateInverseCumulativeDistribution(double x, double p)
{
var n = Normal.WithMeanStdDev(5.0, 2.0);
- AssertHelpers.AlmostEqualRelative(x, n.InverseCumulativeDistribution(f), 14);
- AssertHelpers.AlmostEqualRelative(x, Normal.InvCDF(5.0, 2.0, f), 14);
+ AssertHelpers.AlmostEqualRelative(x, n.InverseCumulativeDistribution(p), 14);
+ AssertHelpers.AlmostEqualRelative(x, Normal.InvCDF(5.0, 2.0, p), 14);
}
///
diff --git a/src/UnitTests/DistributionTests/Continuous/StudentTTests.cs b/src/UnitTests/DistributionTests/Continuous/StudentTTests.cs
index 0f716cb6..8fb54c40 100644
--- a/src/UnitTests/DistributionTests/Continuous/StudentTTests.cs
+++ b/src/UnitTests/DistributionTests/Continuous/StudentTTests.cs
@@ -411,14 +411,6 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
ied.Take(5).ToArray();
}
- ///
- /// Validate cumulative distribution.
- ///
- /// Location value.
- /// Scale value.
- /// Degrees of freedom.
- /// Input X value.
- /// Expected value.
[TestCase(0.0, 1.0, 1.0, 0.0, 0.5)]
[TestCase(0.0, 1.0, 1.0, 1.0, 0.75)]
[TestCase(0.0, 1.0, 1.0, -1.0, 0.25)]
@@ -432,10 +424,31 @@ namespace MathNet.Numerics.UnitTests.DistributionTests.Continuous
[TestCase(0.0, 1.0, Double.PositiveInfinity, 0.0, 0.5)]
[TestCase(0.0, 1.0, Double.PositiveInfinity, 1.0, 0.841344746068543)]
[TestCase(0.0, 1.0, Double.PositiveInfinity, 2.0, 0.977249868051821)]
- public void ValidateCumulativeDistribution(double location, double scale, double dof, double x, double c)
+ public void ValidateCumulativeDistribution(double location, double scale, double dof, double x, double p)
{
- var n = new StudentT(location, scale, dof);
- AssertHelpers.AlmostEqualRelative(c, n.CumulativeDistribution(x), 13);
+ var dist = new StudentT(location, scale, dof);
+ Assert.That(dist.CumulativeDistribution(x), Is.EqualTo(p).Within(1e-13));
+ Assert.That(StudentT.CDF(location, scale, dof, x), Is.EqualTo(p).Within(1e-13));
+ }
+
+ [TestCase(0.0, 1.0, 1.0, 0.0, 0.5)]
+ [TestCase(0.0, 1.0, 1.0, 1.0, 0.75)]
+ [TestCase(0.0, 1.0, 1.0, -1.0, 0.25)]
+ [TestCase(0.0, 1.0, 1.0, 2.0, 0.852416382349567)]
+ [TestCase(0.0, 1.0, 1.0, -2.0, 0.147583617650433)]
+ [TestCase(0.0, 1.0, 2.0, 0.0, 0.5)]
+ [TestCase(0.0, 1.0, 2.0, 1.0, 0.788675134594813)]
+ [TestCase(0.0, 1.0, 2.0, -1.0, 0.211324865405187)]
+ [TestCase(0.0, 1.0, 2.0, 2.0, 0.908248290463863)]
+ [TestCase(0.0, 1.0, 2.0, -2.0, 0.091751709536137)]
+ [TestCase(0.0, 1.0, Double.PositiveInfinity, 0.0, 0.5)]
+ [TestCase(0.0, 1.0, Double.PositiveInfinity, 1.0, 0.841344746068543)]
+ [TestCase(0.0, 1.0, Double.PositiveInfinity, 2.0, 0.977249868051821)]
+ public void ValidateInverseCumulativeDistribution(double location, double scale, double dof, double x, double p)
+ {
+ var dist = new StudentT(location, scale, dof);
+ Assert.That(dist.InverseCumulativeDistribution(p), Is.EqualTo(x).Within(1e-6));
+ Assert.That(StudentT.InvCDF(location, scale, dof, p), Is.EqualTo(x).Within(1e-6));
}
}
}