github.com/gopherd/gonum@v0.0.4/stat/distuv/gamma_test.go (about)

     1  // Copyright ©2016 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  	"fmt"
     9  	"math"
    10  	"sort"
    11  	"testing"
    12  
    13  	"math/rand"
    14  
    15  	"github.com/gopherd/gonum/floats/scalar"
    16  )
    17  
    18  func TestGamma(t *testing.T) {
    19  	t.Parallel()
    20  	// Values a comparison with scipy
    21  	for _, test := range []struct {
    22  		x, alpha, want float64
    23  	}{
    24  		{0.9, 0.1, 0.046986817861555757},
    25  		{0.9, 0.01, 0.0045384353289090401},
    26  		{0.45, 0.01, 0.014137035997241795},
    27  	} {
    28  		pdf := Gamma{Alpha: test.alpha, Beta: 1}.Prob(test.x)
    29  		if !scalar.EqualWithinAbsOrRel(pdf, test.want, 1e-10, 1e-10) {
    30  			t.Errorf("Pdf mismatch. Got %v, want %v", pdf, test.want)
    31  		}
    32  	}
    33  	src := rand.New(rand.NewSource(1))
    34  	for i, g := range []Gamma{
    35  		{Alpha: 0.1, Beta: 0.8, Src: src},
    36  		{Alpha: 0.3, Beta: 0.8, Src: src},
    37  		{Alpha: 0.5, Beta: 0.8, Src: src},
    38  		{Alpha: 0.9, Beta: 6, Src: src},
    39  		{Alpha: 0.9, Beta: 500, Src: src},
    40  		{Alpha: 1, Beta: 1, Src: src},
    41  		{Alpha: 1.6, Beta: 0.4, Src: src},
    42  		{Alpha: 2.6, Beta: 1.5, Src: src},
    43  		{Alpha: 5.6, Beta: 0.5, Src: src},
    44  		{Alpha: 30, Beta: 1.7, Src: src},
    45  		{Alpha: 30.2, Beta: 1.7, Src: src},
    46  	} {
    47  		testGamma(t, g, i)
    48  	}
    49  }
    50  
    51  func testGamma(t *testing.T, f Gamma, i int) {
    52  	// TODO(btracey): Replace this when Gamma implements FullDist.
    53  	const (
    54  		tol  = 1e-2
    55  		n    = 1e5
    56  		bins = 50
    57  	)
    58  	x := make([]float64, n)
    59  	generateSamples(x, f)
    60  	sort.Float64s(x)
    61  
    62  	var quadTol float64
    63  
    64  	if f.Alpha < 0.3 {
    65  		// Gamma PDF has a singularity at 0 for alpha < 1,
    66  		// which gets sharper as alpha -> 0.
    67  		quadTol = 0.2
    68  	} else {
    69  		quadTol = tol
    70  	}
    71  	testRandLogProbContinuous(t, i, 0, x, f, quadTol, bins)
    72  	checkMean(t, i, x, f, tol)
    73  	checkVarAndStd(t, i, x, f, 2e-2)
    74  	checkExKurtosis(t, i, x, f, 2e-1)
    75  	switch {
    76  	case f.Alpha < 0.3:
    77  		quadTol = 0.1
    78  	case f.Alpha < 1:
    79  		quadTol = 1e-3
    80  	default:
    81  		quadTol = 1e-10
    82  	}
    83  	checkProbContinuous(t, i, x, 0, math.Inf(1), f, quadTol)
    84  	checkQuantileCDFSurvival(t, i, x, f, 5e-2)
    85  	if f.Alpha < 1 {
    86  		if !math.IsNaN(f.Mode()) {
    87  			t.Errorf("Expected NaN mode for alpha < 1, got %v", f.Mode())
    88  		}
    89  	} else {
    90  		checkMode(t, i, x, f, 2e-1, 1)
    91  	}
    92  	cdfNegX := f.CDF(-0.0001)
    93  	if cdfNegX != 0 {
    94  		t.Errorf("Expected zero CDF for negative argument, got %v", cdfNegX)
    95  	}
    96  	survNegX := f.Survival(-0.0001)
    97  	if survNegX != 1 {
    98  		t.Errorf("Expected survival function of 1 for negative argument, got %v", survNegX)
    99  	}
   100  	if f.NumParameters() != 2 {
   101  		t.Errorf("Mismatch in NumParameters: got %v, want 2", f.NumParameters())
   102  	}
   103  	lPdf := f.LogProb(-0.0001)
   104  	if !math.IsInf(lPdf, -1) {
   105  		t.Errorf("Expected log(CDF) to be -Infinity for negative argument, got %v", lPdf)
   106  	}
   107  }
   108  
   109  func TestGammaPanics(t *testing.T) {
   110  	t.Parallel()
   111  	g := Gamma{1, 0, nil}
   112  	if !panics(func() { g.Rand() }) {
   113  		t.Errorf("Expected Rand panic for Beta <= 0")
   114  	}
   115  	g = Gamma{0, 1, nil}
   116  	if !panics(func() { g.Rand() }) {
   117  		t.Errorf("Expected Rand panic for Alpha <= 0")
   118  	}
   119  }
   120  
   121  func BenchmarkGammaRand(b *testing.B) {
   122  	src := rand.New(rand.NewSource(1))
   123  	for i, g := range []Gamma{
   124  		{Alpha: 0.1, Beta: 0.8, Src: src},
   125  		{Alpha: 0.5, Beta: 0.8, Src: src},
   126  		{Alpha: 1, Beta: 1, Src: src},
   127  		{Alpha: 1.6, Beta: 0.4, Src: src},
   128  		{Alpha: 30, Beta: 1.7, Src: src},
   129  	} {
   130  		b.Run(fmt.Sprintf("case %d", i), func(b *testing.B) {
   131  			for i := 0; i < b.N; i++ {
   132  				g.Rand()
   133  			}
   134  		})
   135  	}
   136  }