Fix #121 and add dependency on IrrationalConstants (extension of #118) (

#122)

Co-authored-by: Tamas K. Papp <tkpapp@gmail.com>
Co-authored-by: Tor Erlend Fjelde <tor.erlend95@gmail.com>

Project.toml

@@ -1,14 +1,18 @@
 name = "StatsFuns"
 uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c"
-version = "0.9.8"
+version = "0.9.9"
 
 [deps]
+IrrationalConstants = "92d709cd-6900-40b7-9082-c6be49f344b6"
 LogExpFunctions = "2ab3a3ac-af41-5b50-aa03-7779005ae688"
+Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 Rmath = "79098fc4-a85e-5d69-aa6a-4863f24498fa"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
 
 [compat]
-LogExpFunctions = "0.2.1"
+IrrationalConstants = "0.1"
+LogExpFunctions = "0.3"
+Reexport = "1"
 Rmath = "0.4, 0.5, 0.6, 0.7"
 SpecialFunctions = "0.8, 0.9, 0.10, 1.0"
 julia = "1"

src/StatsFuns.jl

@@ -2,11 +2,12 @@ __precompile__(true)
 
 module StatsFuns
 
-import Base: Math.@horner, @irrational
+using Base: Math.@horner
+using Reexport
 using SpecialFunctions
 
-export
- # constants
+# reexports
+@reexport using IrrationalConstants:
  twoπ, # 2π
  fourπ, # 4π
  halfπ, # π / 2
@@ -24,7 +25,33 @@ export
  sqrthalfπ, # √(π / 2)
  invsqrt2, # 1 / √2
  invsqrt2π, # 1 / √2π
+ loghalf, # log(1 / 2)
+ logtwo, # log(2)
+ logπ, # log(π)
+ log2π, # log(2π)
+ log4π # log(4π)
 
+@reexport using LogExpFunctions:
+ xlogx, # x * log(x) for x > 0, or 0 when x == 0
+ xlogy, # x * log(y) for x > 0, or 0 when x == 0
+ logistic, # 1 / (1 + exp(-x))
+ logit, # log(x / (1 - x))
+ log1psq, # log(1 + x^2)
+ log1pexp, # log(1 + exp(x))
+ log1mexp, # log(1 - exp(x))
+ log2mexp, # log(2 - exp(x))
+ logexpm1, # log(exp(x) - 1)
+ softplus, # alias of log1pexp
+ invsoftplus, # alias of logexpm1
+ log1pmx, # log(1 + x) - x
+ logmxp1, # log(x) - x + 1
+ logaddexp, # log(exp(x) + exp(y))
+ logsubexp, # log(abs(e^x - e^y))
+ logsumexp, # log(sum(exp(x)))
+ softmax, # exp(x_i) / sum(exp(x)), for i
+ softmax! # inplace softmax
+
+export
  # distrs/beta
  betapdf, # pdf of beta distribution
  betalogpdf, # logpdf of beta distribution
@@ -208,36 +235,7 @@ export
  logmvbeta, # logarithm of multivariate beta function
  lstirling_asym
 
-# reexports
-using LogExpFunctions
-
-export
- loghalf, # log(1 / 2)
- logtwo, # log(2)
- logπ, # log(π)
- log2π, # log(2π)
- log4π, # log(4π)
- xlogx, # x * log(x) for x > 0, or 0 when x == 0
- xlogy, # x * log(y) for x > 0, or 0 when x == 0
- logistic, # 1 / (1 + exp(-x))
- logit, # log(x / (1 - x))
- log1psq, # log(1 + x^2)
- log1pexp, # log(1 + exp(x))
- log1mexp, # log(1 - exp(x))
- log2mexp, # log(2 - exp(x))
- logexpm1, # log(exp(x) - 1)
- softplus, # alias of log1pexp
- invsoftplus, # alias of logexpm1
- log1pmx, # log(1 + x) - x
- logmxp1, # log(x) - x + 1
- logaddexp, # log(exp(x) + exp(y))
- logsubexp, # log(abs(e^x - e^y))
- logsumexp, # log(sum(exp(x)))
- softmax, # exp(x_i) / sum(exp(x)), for i
- softmax! # inplace softmax
-
 ## source files
-include("constants.jl")
 include("misc.jl")
 include("rmath.jl")
 

src/misc.jl

@@ -7,7 +7,7 @@ Return the logarithm of [multivariate gamma function](https://en.wikipedia.org/w
 """
 function logmvgamma(p::Int, a::Real)
  # NOTE: one(a) factors are here to prevent unnecessary promotion of Float32
- res = p * (p - 1) * log(pi * one(a)) / 4
+ res = p * (p - 1) * (logπ * one(a)) / 4
  for ii in 1:p
  res += loggamma(a + (1 - ii) * one(a) / 2)
  end
@@ -19,7 +19,8 @@ end
 
 Return the logarithm of the multivariate beta function ([DLMF 35.3.7](https://dlmf.nist.gov/35.3#E7)).
 """
-logmvbeta(p::Int, a::Real, b::Real) = logmvgamma(p, a) + logmvgamma(p, b) - logmvgamma(p, a + b)
+logmvbeta(p::Int, a::T, b::T) where {T<:Real} = logmvgamma(p, a) + logmvgamma(p, b) - logmvgamma(p, a + b)
+logmvbeta(p::Int, a::Real, b::Real) = logmvbeta(p, promote(a, b)...)
 
 """
  lstirling_asym(x)