Remove code for unsupported Julia versions (#155)

* Remove code for unsupported Julia versions

* Partially revert change

* Slightly improve `_expm1(::Float16)`

* Another fix

* Stay in log space

* Another typo...

* Add `_log2mexp`

* Just fall back to Float32

* Always use `LogExpFunctions.log2mexp`

Project.toml

@@ -25,7 +25,7 @@ ChainRulesCore = "1"
 HypergeometricFunctions = "0.3.10"
 InverseFunctions = "0.1"
 IrrationalConstants = "0.1, 0.2"
-LogExpFunctions = "0.3.2"
+LogExpFunctions = "0.3.26"
 Reexport = "1"
 Rmath = "0.6.1, 0.7"
 SpecialFunctions = "1.8.4, 2.1.4"

src/distrs/tdist.jl

@@ -48,25 +48,16 @@ end
 tdistinvcdf(ν::Real, p::Real) = tdistinvcdf(map(float, promote(ν, p))...)
 
 tdistinvccdf(ν::Real, p::Real) = -tdistinvcdf(ν, p)
-
-if VERSION < v"1.7.0-DEV.1172"
-    function _expm1(x::Float16)
-        if -0.2 < x < 0.1
-            return @evalpoly(x, Float16(0), Float16(1), Float16(1/2), Float16(1/6), Float16(1/24), Float16(1/120))
-        else
-            return exp(x) - 1
-        end
-    end
-end
-_expm1(x::Number) = expm1(x)
-
 function tdistinvlogcdf(ν::T, logp::T) where T<:Real
     if isinf(ν)
         return norminvlogcdf(logp)
-    elseif logp < -log(2)
-        return -sqrt(fdistinvlogccdf(one(ν), ν, logp + logtwo))
     else
-        return sqrt(fdistinvccdf(one(ν), ν, -2*_expm1(logp)))
+        logq = logp + logtwo
+        if logq < 0
+            return -sqrt(fdistinvlogccdf(one(ν), ν, logq))
+        else
+            return sqrt(fdistinvlogccdf(one(ν), ν, log2mexp(logq)))
+        end
     end
 end
 tdistinvlogcdf(ν::Real, logp::Real) = tdistinvlogcdf(map(float, promote(ν, logp))...)

test/misc.jl

@@ -43,12 +43,7 @@ end
             #  B₁(a, b) = B(a, b)
             a = rand(eltya)
             b = rand(eltyb)
-            # Older SpecialFunctions + Julia versions require slightly larger tolerance
-            if VERSION < v"1.3" && promote_type(eltya, eltyb) === Float64
-                @test logmvbeta(1, a, b) ≈ logbeta(a, b) rtol = 10 * sqrt(eps(Float64))
-            else
-                @test logmvbeta(1, a, b) ≈ logbeta(a, b)
-            end
+            @test logmvbeta(1, a, b) ≈ logbeta(a, b)
         end
     end