stat/distuv: Add quantile and survival functions to Gamma, ChiSquared…

…, and F

distuv/chisquared.go

@@ -72,11 +72,27 @@ func (c ChiSquared) Rand() float64 {
 	return Gamma{c.K / 2, 0.5, c.Src}.Rand()
 }
 
+// Quantile returns the inverse of the cumulative distribution function.
+func (c ChiSquared) Quantile(p float64) float64 {
+	if p < 0 || p > 1 {
+		panic(badPercentile)
+	}
+	return mathext.GammaIncInv(0.5*c.K, p) * 2
+}
+
 // StdDev returns the standard deviation of the probability distribution.
 func (c ChiSquared) StdDev() float64 {
 	return math.Sqrt(c.Variance())
 }
 
+// Survival returns the survival function (complementary CDF) at x.
+func (c ChiSquared) Survival(x float64) float64 {
+	if x < 0 {
+		return 1
+	}
+	return mathext.GammaIncComp(0.5*c.K, 0.5*x)
+}
+
 // Variance returns the variance of the probability distribution.
 func (c ChiSquared) Variance() float64 {
 	return 2 * c.K

distuv/chisquared_test.go

@@ -74,4 +74,5 @@ func testChiSquared(t *testing.T, c ChiSquared, i int) {
 	checkVarAndStd(t, i, x, c, tol)
 	checkExKurtosis(t, i, x, c, 5e-2)
 	checkProbContinuous(t, i, x, c, 1e-3)
+	checkQuantileCDFSurvival(t, i, x, c, 1e-3)
 }

distuv/f.go

@@ -76,6 +76,15 @@ func (f F) Prob(x float64) float64 {
 	return math.Exp(f.LogProb(x))
 }
 
+// Quantile returns the inverse of the cumulative distribution function.
+func (f F) Quantile(p float64) float64 {
+	if p < 0 || p > 1 {
+		panic(badPercentile)
+	}
+	y := mathext.InvRegIncBeta(0.5*f.D1, 0.5*f.D2, p)
+	return f.D2 * y / (f.D1 * (1 - y))
+}
+
 // Rand returns a random sample drawn from the distribution.
 func (f F) Rand() float64 {
 	u1 := ChiSquared{f.D1, f.Source}.Rand()
@@ -105,6 +114,11 @@ func (f F) StdDev() float64 {
 	return math.Sqrt(f.Variance())
 }
 
+// Survival returns the survival function (complementary CDF) at x.
+func (f F) Survival(x float64) float64 {
+	return 1 - f.CDF(x)
+}
+
 // Variance returns the variance of the probability distribution.
 //
 // Variance returns NaN if the D2 parameter is less than or equal to 4.

distuv/f_test.go

@@ -85,4 +85,5 @@ func testF(t *testing.T, f F, i int) {
 	checkVarAndStd(t, i, x, f, tol)
 	checkExKurtosis(t, i, x, f, 5e-2)
 	checkSkewness(t, i, x, f, tol)
+	checkQuantileCDFSurvival(t, i, x, f, 1e-3)
 }

distuv/gamma.go

@@ -7,6 +7,8 @@ package distuv
 import (
 	"math"
 	"math/rand"
+
+	"github.com/gonum/mathext"
 )
 
 // Gamma implements the Gamma distribution, a two-parameter continuous distribution
@@ -27,6 +29,14 @@ type Gamma struct {
 	Source *rand.Rand
 }
 
+// CDF computes the value of the cumulative distribution function at x.
+func (g Gamma) CDF(x float64) float64 {
+	if x < 0 {
+		return 0
+	}
+	return mathext.GammaInc(g.Alpha, g.Beta*x)
+}
+
 // ExKurtosis returns the excess kurtosis of the distribution.
 func (g Gamma) ExKurtosis() float64 {
 	return 6 / g.Alpha
@@ -70,6 +80,14 @@ func (g Gamma) Prob(x float64) float64 {
 	return math.Exp(g.LogProb(x))
 }
 
+// Quantile returns the inverse of the cumulative distribution function.
+func (g Gamma) Quantile(p float64) float64 {
+	if p < 0 || p > 1 {
+		panic(badPercentile)
+	}
+	return mathext.GammaIncInv(g.Alpha, p) / g.Beta
+}
+
 // Rand returns a random sample drawn from the distribution.
 //
 // Rand panics if either alpha or beta is <= 0.
@@ -207,6 +225,14 @@ func (g Gamma) Rand() float64 {
 	panic("unreachable")
 }
 
+// Survival returns the survival function (complementary CDF) at x.
+func (g Gamma) Survival(x float64) float64 {
+	if x < 0 {
+		return 1
+	}
+	return mathext.GammaIncComp(g.Alpha, g.Beta*x)
+}
+
 // StdDev returns the standard deviation of the probability distribution.
 func (g Gamma) StdDev() float64 {
 	return math.Sqrt(g.Variance())

distuv/gamma_test.go

@@ -59,4 +59,5 @@ func testGamma(t *testing.T, f Gamma, i int) {
 	checkVarAndStd(t, i, x, f, 2e-2)
 	checkExKurtosis(t, i, x, f, 5e-2)
 	checkProbContinuous(t, i, x, f, 1e-3)
+	checkQuantileCDFSurvival(t, i, x, f, 1e-2)
 }