From 279075ae2d80e2c98d7728b001b5a713a73e0bb8 Mon Sep 17 00:00:00 2001
From: Andreas Noack <andreas@noack.dk>
Date: Tue, 13 Feb 2024 21:47:38 +0100
Subject: [PATCH] =?UTF-8?q?Handle=20tiny=20=CE=B8=20in=20gammalogpdf?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

---
 src/distrs/gamma.jl | 7 ++++++-
 test/rmath.jl       | 1 +
 2 files changed, 7 insertions(+), 1 deletion(-)

diff --git a/src/distrs/gamma.jl b/src/distrs/gamma.jl
index 2ed1bde..abc0f29 100644
--- a/src/distrs/gamma.jl
+++ b/src/distrs/gamma.jl
@@ -22,7 +22,12 @@ gammalogpdf(k::Real, θ::Real, x::Real) = gammalogpdf(promote(k, θ, x)...)
 function gammalogpdf(k::T, θ::T, x::T) where {T<:Real}
     # we ensure that `log(x)` does not error if `x < 0`
     xθ = max(x, 0) / θ
-    val = -loggamma(k) + xlogy(k - 1, xθ) - log(θ) - xθ
+    val = -loggamma(k) - log(θ) - xθ
+    # xlogy(k - 1, xθ) - xθ -> -∞ for xθ -> ∞ so we only add the first term
+    # when it's safe
+    if isfinite(xθ)
+        val += xlogy(k - 1, xθ)
+    end
     return x < 0 ? oftype(val, -Inf) : val
 end
 
diff --git a/test/rmath.jl b/test/rmath.jl
index 4582ffa..0276353 100644
--- a/test/rmath.jl
+++ b/test/rmath.jl
@@ -264,6 +264,7 @@ end
         ((Float16(1), Float16(1)), (Float16(0):Float16(0.05):Float16(12))),
         ((1f0, 1f0), (Float16(0):Float16(0.05):Float16(12))),
         ((2, 3), (0//1:12//1)),
+        ((3.0, 1e-310), (0.0:0.05:12.0))
     ])
 
     rmathcomp_tests("hyper", [