gonum.org/v1/gonum@v0.14.0/stat/distuv/f.go (about)

     1  // Copyright ©2017 The Gonum Authors. All rights reserved.
     2  // Use of this source code is governed by a BSD-style
     3  // license that can be found in the LICENSE file.
     4  
     5  package distuv
     6  
     7  import (
     8  	"math"
     9  
    10  	"golang.org/x/exp/rand"
    11  
    12  	"gonum.org/v1/gonum/mathext"
    13  )
    14  
    15  // F implements the F-distribution, a two-parameter continuous distribution
    16  // with support over the positive real numbers.
    17  //
    18  // The F-distribution has density function
    19  //
    20  //	sqrt(((d1*x)^d1) * d2^d2 / ((d1*x+d2)^(d1+d2))) / (x * B(d1/2,d2/2))
    21  //
    22  // where B is the beta function.
    23  //
    24  // For more information, see https://en.wikipedia.org/wiki/F-distribution
    25  type F struct {
    26  	D1  float64 // Degrees of freedom for the numerator
    27  	D2  float64 // Degrees of freedom for the denominator
    28  	Src rand.Source
    29  }
    30  
    31  // CDF computes the value of the cumulative density function at x.
    32  func (f F) CDF(x float64) float64 {
    33  	return mathext.RegIncBeta(f.D1/2, f.D2/2, f.D1*x/(f.D1*x+f.D2))
    34  }
    35  
    36  // ExKurtosis returns the excess kurtosis of the distribution.
    37  //
    38  // ExKurtosis returns NaN if the D2 parameter is less or equal to 8.
    39  func (f F) ExKurtosis() float64 {
    40  	if f.D2 <= 8 {
    41  		return math.NaN()
    42  	}
    43  	return (12 / (f.D2 - 6)) * ((5*f.D2-22)/(f.D2-8) + ((f.D2-4)/f.D1)*((f.D2-2)/(f.D2-8))*((f.D2-2)/(f.D1+f.D2-2)))
    44  }
    45  
    46  // LogProb computes the natural logarithm of the value of the probability
    47  // density function at x.
    48  func (f F) LogProb(x float64) float64 {
    49  	return 0.5*(f.D1*math.Log(f.D1*x)+f.D2*math.Log(f.D2)-(f.D1+f.D2)*math.Log(f.D1*x+f.D2)) - math.Log(x) - mathext.Lbeta(f.D1/2, f.D2/2)
    50  }
    51  
    52  // Mean returns the mean of the probability distribution.
    53  //
    54  // Mean returns NaN if the D2 parameter is less than or equal to 2.
    55  func (f F) Mean() float64 {
    56  	if f.D2 <= 2 {
    57  		return math.NaN()
    58  	}
    59  	return f.D2 / (f.D2 - 2)
    60  }
    61  
    62  // Mode returns the mode of the distribution.
    63  //
    64  // Mode returns NaN if the D1 parameter is less than or equal to 2.
    65  func (f F) Mode() float64 {
    66  	if f.D1 <= 2 {
    67  		return math.NaN()
    68  	}
    69  	return ((f.D1 - 2) / f.D1) * (f.D2 / (f.D2 + 2))
    70  }
    71  
    72  // NumParameters returns the number of parameters in the distribution.
    73  func (f F) NumParameters() int {
    74  	return 2
    75  }
    76  
    77  // Prob computes the value of the probability density function at x.
    78  func (f F) Prob(x float64) float64 {
    79  	return math.Exp(f.LogProb(x))
    80  }
    81  
    82  // Quantile returns the inverse of the cumulative distribution function.
    83  func (f F) Quantile(p float64) float64 {
    84  	if p < 0 || p > 1 {
    85  		panic(badPercentile)
    86  	}
    87  	y := mathext.InvRegIncBeta(0.5*f.D1, 0.5*f.D2, p)
    88  	return f.D2 * y / (f.D1 * (1 - y))
    89  }
    90  
    91  // Rand returns a random sample drawn from the distribution.
    92  func (f F) Rand() float64 {
    93  	u1 := ChiSquared{f.D1, f.Src}.Rand()
    94  	u2 := ChiSquared{f.D2, f.Src}.Rand()
    95  	return (u1 / f.D1) / (u2 / f.D2)
    96  }
    97  
    98  // Skewness returns the skewness of the distribution.
    99  //
   100  // Skewness returns NaN if the D2 parameter is less than or equal to 6.
   101  func (f F) Skewness() float64 {
   102  	if f.D2 <= 6 {
   103  		return math.NaN()
   104  	}
   105  	num := (2*f.D1 + f.D2 - 2) * math.Sqrt(8*(f.D2-4))
   106  	den := (f.D2 - 6) * math.Sqrt(f.D1*(f.D1+f.D2-2))
   107  	return num / den
   108  }
   109  
   110  // StdDev returns the standard deviation of the probability distribution.
   111  //
   112  // StdDev returns NaN if the D2 parameter is less than or equal to 4.
   113  func (f F) StdDev() float64 {
   114  	if f.D2 <= 4 {
   115  		return math.NaN()
   116  	}
   117  	return math.Sqrt(f.Variance())
   118  }
   119  
   120  // Survival returns the survival function (complementary CDF) at x.
   121  func (f F) Survival(x float64) float64 {
   122  	return 1 - f.CDF(x)
   123  }
   124  
   125  // Variance returns the variance of the probability distribution.
   126  //
   127  // Variance returns NaN if the D2 parameter is less than or equal to 4.
   128  func (f F) Variance() float64 {
   129  	if f.D2 <= 4 {
   130  		return math.NaN()
   131  	}
   132  	num := 2 * f.D2 * f.D2 * (f.D1 + f.D2 - 2)
   133  	den := f.D1 * (f.D2 - 2) * (f.D2 - 2) * (f.D2 - 4)
   134  	return num / den
   135  }