Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Provides same interface to studentized range distribution functions ([p/q]tukey) #66

Merged
merged 13 commits into from
Feb 4, 2019
Merged

Provides same interface to studentized range distribution functions ([p/q]tukey) #66

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
13 commits
Select commit Hold shift + click to select a range
29b4076
Added studentized range distribution (srdist)
BioTurboNick Jan 25, 2019
e61cc29
Fixed interface to libRmath
BioTurboNick Jan 25, 2019
7ff8304
Partial tests, but need to better understand what they're doing
BioTurboNick Jan 25, 2019
85e1718
Wrapped in TestSet, documented case where functions diverge
BioTurboNick Jan 25, 2019
1e8e3af
Removed stray bracket
BioTurboNick Jan 26, 2019
694dbde
Corrected pdf
BioTurboNick Jan 27, 2019
2487075
Simplified pdf code
BioTurboNick Jan 27, 2019
5e0d949
Switched srdistpdf internal functions to use normpdf/cdf
BioTurboNick Jan 29, 2019
dbc1d63
Removed ForwardDiff
BioTurboNick Jan 29, 2019
34a4b92
Wrapped generic tests in TestSet too
BioTurboNick Jan 29, 2019
087c0aa
Removed pdf and associated tests
BioTurboNick Jan 29, 2019
bc37581
Removed QuadGK requirement
BioTurboNick Jan 29, 2019
259bcda
Consolidated tukey cases with others
BioTurboNick Jan 29, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
11 changes: 11 additions & 0 deletions src/StatsFuns.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -217,6 +217,16 @@ export
tdistinvlogcdf, # inverse-logcdf of student's t distribution
tdistinvlogccdf, # inverse-logccdf of student's t distribution

# distrs/srdist
srdistcdf, # cdf of studentized range distribution
srdistccdf, # ccdf of studentized range distribution
srdistlogcdf, # logcdf of studentized range distribution
srdistlogccdf, # logccdf of studentized range distribution
srdistinvcdf, # inverse-cdf of studentized range distribution
srdistinvccdf, # inverse-ccdf of studentized range distribution
srdistinvlogcdf, # inverse-logcdf of studentized range distribution
srdistinvlogccdf, # inverse-logccdf of studentized range distribution

# misc
logmvgamma, # logarithm of multivariate gamma function
lstirling_asym
Expand Down Expand Up @@ -244,5 +254,6 @@ include(joinpath("distrs", "norm.jl"))
include(joinpath("distrs", "ntdist.jl"))
include(joinpath("distrs", "pois.jl"))
include(joinpath("distrs", "tdist.jl"))
include(joinpath("distrs", "srdist.jl"))

end # module
4 changes: 2 additions & 2 deletions src/distrs/fdist.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ import .RFunctions:
fdistinvlogccdf

# pdf for numbers with generic types
fdistpdf(d1::Real, d2::Real, x::Number) = sqrt((d1 * x)^d1 * d2^d2 / (d1 * x + d2)^(d1 + d2)) / (x * beta(d1 / 2, d2 / 2))
fdistpdf(ν1::Real, ν2::Real, x::Number) = sqrt((ν1 * x)^ν1 * ν2^ν2 / (ν1 * x + ν2)^(ν1 + ν2)) / (x * beta(ν1 / 2, ν2 / 2))

# logpdf for numbers with generic types
fdistlogpdf(d1::Real, d2::Real, x::Number) = (d1 * log(d1 * x) + d2 * log(d2) - (d1 + d2) * log(d1 * x + d2)) / 2 - log(x) - lbeta(d1 / 2, d2 / 2)
fdistlogpdf(ν1::Real, ν2::Real, x::Number) = (ν1 * log(ν1 * x) + ν2 * log(ν2) - (ν1 + ν2) * log(ν1 * x + ν2)) / 2 - log(x) - lbeta(ν1 / 2, ν2 / 2)
11 changes: 11 additions & 0 deletions src/distrs/srdist.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,11 @@
# functions related to studentized range distribution

import .RFunctions:
srdistcdf,
srdistccdf,
srdistlogcdf,
srdistlogccdf,
srdistinvcdf,
srdistinvccdf,
srdistinvlogcdf,
srdistinvlogccdf
39 changes: 28 additions & 11 deletions src/rmath.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -34,29 +34,45 @@ function _import_rmath(rname::Symbol, jname::Symbol, pargs)
invlogcdf = Symbol(jname, "invlogcdf")
invlogccdf = Symbol(jname, "invlogccdf")

is_tukey = rname == :tukey

rand = Symbol(jname, "rand")
has_rand = true
if rname == :nbeta || rname == :nf || rname == :nt
if rname == :nbeta || rname == :nf || rname == :nt || is_tukey
has_rand = false
end

# arguments & argument types
_pts = fill(:Cdouble, length(pargs))

dtypes = Expr(:tuple, :Cdouble, _pts..., :Cint)
ptypes = Expr(:tuple, :Cdouble, _pts..., :Cint, :Cint)
qtypes = Expr(:tuple, :Cdouble, _pts..., :Cint, :Cint)
rtypes = Expr(:tuple, _pts...)
if is_tukey
dtypes = Expr(:tuple, :Cdouble, :Cdouble, _pts..., :Cint)
ptypes = Expr(:tuple, :Cdouble, :Cdouble, _pts..., :Cint, :Cint)
qtypes = Expr(:tuple, :Cdouble, :Cdouble, _pts..., :Cint, :Cint)
rtypes = Expr(:tuple, :Cdouble, _pts...)
else
dtypes = Expr(:tuple, :Cdouble, _pts..., :Cint)
ptypes = Expr(:tuple, :Cdouble, _pts..., :Cint, :Cint)
qtypes = Expr(:tuple, :Cdouble, _pts..., :Cint, :Cint)
rtypes = Expr(:tuple, _pts...)
end

pdecls = [Expr(:(::), ps, :(Union{Float64,Int})) for ps in pargs] # [:(p1::Union{Float64, Int}), :(p2::Union{...}), ...]

# function implementation
if is_tukey
# ptukey and qtukey have an extra literal 1 argument
pargs = (1, pargs...)
end

# Function implementation
quote
$pdf($(pdecls...), x::Union{Float64,Int}) =
ccall(($dfun, libRmath), Float64, $dtypes, x, $(pargs...), 0)
if $(!is_tukey)
$pdf($(pdecls...), x::Union{Float64,Int}) =
ccall(($dfun, libRmath), Float64, $dtypes, x, $(pargs...), 0)

$logpdf($(pdecls...), x::Union{Float64,Int}) =
ccall(($dfun, libRmath), Float64, $dtypes, x, $(pargs...), 1)
$logpdf($(pdecls...), x::Union{Float64,Int}) =
ccall(($dfun, libRmath), Float64, $dtypes, x, $(pargs...), 1)
end

$cdf($(pdecls...), x::Union{Float64,Int}) =
ccall(($pfun, libRmath), Float64, $ptypes, x, $(pargs...), 1, 0)
Expand Down Expand Up @@ -99,13 +115,14 @@ end
@import_rmath beta beta α β
@import_rmath binom binom n p
@import_rmath chisq chisq k
@import_rmath f fdist d1 d2
@import_rmath f fdist ν1 ν2
@import_rmath gamma gamma α β
@import_rmath hyper hyper ms mf n
@import_rmath norm norm μ σ
@import_rmath nbinom nbinom r p
@import_rmath pois pois λ
@import_rmath t tdist k
@import_rmath tukey srdist k ν

@import_rmath nbeta nbeta α β λ
@import_rmath nchisq nchisq k λ
Expand Down
174 changes: 92 additions & 82 deletions test/generic.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -43,85 +43,95 @@ end

### Test cases

genericcomp_tests("beta", [
((1.0, 1.0), 0.01:0.01:0.99),
((2.0, 3.0), 0.01:0.01:0.99),
((10.0, 2.0), 0.01:0.01:0.99),
])

genericcomp_tests("binom", [
((1, 0.5), 0.0:1.0),
((1, 0.7), 0.0:1.0),
((8, 0.6), 0.0:8.0),
((20, 0.1), 0.0:20.0),
((20, 0.9), 0.0:20.0),
((20, 0.9), 0:20),
])

genericcomp_tests("chisq", [
((1,), 0.0:0.1:8.0),
((4,), 0.0:0.1:8.0),
((9,), 0.0:0.1:8.0),
])

genericcomp_tests("fdist", [
((1, 1), (0.0:0.1:5.0)),
((2, 1), (0.0:0.1:5.0)),
((5, 2), (0.0:0.1:5.0)),
((10, 1), (0.0:0.1:5.0)),
((10, 3), (0.0:0.1:5.0)),
])

genericcomp_tests("gamma", [
((1.0, 1.0), (0.05:0.05:12.0)),
((0.5, 1.0), (0.05:0.05:12.0)),
((3.0, 1.0), (0.05:0.05:12.0)),
((9.0, 1.0), (0.05:0.05:12.0)),
((2.0, 3.0), (0.05:0.05:12.0)),
])

# genericcomp_tests("nbeta", [
# ((1.0, 1.0, 0.0), 0.01:0.01:0.99),
# ((2.0, 3.0, 0.0), 0.01:0.01:0.99),
# ((1.0, 1.0, 2.0), 0.01:0.01:0.99),
# ((3.0, 4.0, 2.0), 0.01:0.01:0.99),
# ])

# genericcomp_tests("nchisq", [
# ((2, 1), 0.0:0.2:8.0),
# ((2, 3), 0.0:0.2:8.0),
# ((4, 1), 0.0:0.2:8.0),
# ((4, 3), 0.0:0.2:8.0),
# ])

# genericcomp_tests("nfdist", [
# ((1.0, 1.0, 0.0), 0.1:0.1:10.0),
# ((1.0, 1.0, 2.0), 0.1:0.1:10.0),
# ((2.0, 3.0, 1.0), 0.1:0.1:10.0),
# ])

genericcomp_tests("norm", [
((0.0, 1.0), -6.0:0.01:6.0),
((2.0, 1.0), -3.0:0.01:7.0),
((0.0, 0.5), -3.0:0.01:3.0),
((0, 1), -3.0:3.0),
((0, 1), -3.0:0.01:3.0)
])

# genericcomp_tests("ntdist", [
# ((0, 1), -4.0:0.1:10.0),
# ((0, 4), -4.0:0.1:10.0),
# ((2, 1), -4.0:0.1:10.0),
# ((2, 4), -4.0:0.1:10.0),
# ])

genericcomp_tests("pois", [
((1.0,), 0:30),
((10.0,), 0:42)
])

genericcomp_tests("tdist", [
((1,), -5.0:0.1:5.0),
((2,), -5.0:0.1:5.0),
((5,), -5.0:0.1:5.0),
])
@testset "Generic" begin
genericcomp_tests("beta", [
((1.0, 1.0), 0.01:0.01:0.99),
((2.0, 3.0), 0.01:0.01:0.99),
((10.0, 2.0), 0.01:0.01:0.99),
])

genericcomp_tests("binom", [
((1, 0.5), 0.0:1.0),
((1, 0.7), 0.0:1.0),
((8, 0.6), 0.0:8.0),
((20, 0.1), 0.0:20.0),
((20, 0.9), 0.0:20.0),
((20, 0.9), 0:20),
])

genericcomp_tests("chisq", [
((1,), 0.0:0.1:8.0),
((4,), 0.0:0.1:8.0),
((9,), 0.0:0.1:8.0),
])

genericcomp_tests("fdist", [
((1, 1), (0.0:0.1:5.0)),
((2, 1), (0.0:0.1:5.0)),
((5, 2), (0.0:0.1:5.0)),
((10, 1), (0.0:0.1:5.0)),
((10, 3), (0.0:0.1:5.0)),
])

genericcomp_tests("gamma", [
((1.0, 1.0), (0.05:0.05:12.0)),
((0.5, 1.0), (0.05:0.05:12.0)),
((3.0, 1.0), (0.05:0.05:12.0)),
((9.0, 1.0), (0.05:0.05:12.0)),
((2.0, 3.0), (0.05:0.05:12.0)),
])

# genericcomp_tests("nbeta", [
# ((1.0, 1.0, 0.0), 0.01:0.01:0.99),
# ((2.0, 3.0, 0.0), 0.01:0.01:0.99),
# ((1.0, 1.0, 2.0), 0.01:0.01:0.99),
# ((3.0, 4.0, 2.0), 0.01:0.01:0.99),
# ])

# genericcomp_tests("nchisq", [
# ((2, 1), 0.0:0.2:8.0),
# ((2, 3), 0.0:0.2:8.0),
# ((4, 1), 0.0:0.2:8.0),
# ((4, 3), 0.0:0.2:8.0),
# ])

# genericcomp_tests("nfdist", [
# ((1.0, 1.0, 0.0), 0.1:0.1:10.0),
# ((1.0, 1.0, 2.0), 0.1:0.1:10.0),
# ((2.0, 3.0, 1.0), 0.1:0.1:10.0),
# ])

genericcomp_tests("norm", [
((0.0, 1.0), -6.0:0.01:6.0),
((2.0, 1.0), -3.0:0.01:7.0),
((0.0, 0.5), -3.0:0.01:3.0),
((0, 1), -3.0:3.0),
((0, 1), -3.0:0.01:3.0)
])

# genericcomp_tests("ntdist", [
# ((0, 1), -4.0:0.1:10.0),
# ((0, 4), -4.0:0.1:10.0),
# ((2, 1), -4.0:0.1:10.0),
# ((2, 4), -4.0:0.1:10.0),
# ])

genericcomp_tests("pois", [
((1.0,), 0:30),
((10.0,), 0:42)
])

genericcomp_tests("tdist", [
((1,), -5.0:0.1:5.0),
((2,), -5.0:0.1:5.0),
((5,), -5.0:0.1:5.0),
])

#genericcomp_tests("srdist", [
# ((1,2), (0.0:0.1:5.0)),
# ((2,2), (0.0:0.1:5.0)),
# ((5,3), (0.0:0.1:5.0)),
# ((10,2), (0.0:0.1:5.0)),
# ((10,5), (0.0:0.1:5.0))
#])
end
Loading