Merge pull request #99 from juliangehring/julia-v0.7

Support for julia v0.7

Author: juliangehring

.appveyor.yml

@@ -1,12 +1,9 @@
 environment:
   matrix:
-  - JULIA_URL: "https://julialang-s3.julialang.org/bin/winnt/x64/0.6/julia-0.6-latest-win64.exe"
-  - JULIA_URL: "https://julialangnightlies-s3.julialang.org/bin/winnt/x64/julia-latest-win64.exe"
+  - JULIA_URL: "https://julialang-s3.julialang.org/bin/winnt/x64/0.7/julia-0.7-latest-win64.exe"
 
 matrix:
   fast_finish: true
-  allow_failures:
-    - JULIA_URL: "https://julialangnightlies-s3.julialang.org/bin/winnt/x64/julia-latest-win64.exe"
 
 notifications:
   - provider: Email

.travis.yml

@@ -4,13 +4,10 @@ os:
   - linux
 
 julia:
-  - 0.6
-  - nightly
+  - 0.7
 
 matrix:
   fast_finish: true
-  allow_failures:
-    - julia: nightly
 
 notifications:
   email:
@@ -23,5 +20,4 @@ script:
 
 after_success:
   - julia -e 'cd(Pkg.dir("MultipleTesting")); Pkg.add("Coverage"); using Coverage; Codecov.submit(Codecov.process_folder(), flags="julia_$(VERSION.major)_$(VERSION.minor)")'
-  - julia -e 'Pkg.add("Documenter")'
-  - julia -e 'cd(Pkg.dir("MultipleTesting")); include(joinpath("docs", "make.jl"))'
+  - julia -e 'cd(Pkg.dir("MultipleTesting")); Pkg.add("Documenter"); include(joinpath("docs", "make.jl"))'

REQUIRE

@@ -1,4 +1,4 @@
-julia 0.6
-StatsBase 0.15.0
-Distributions 0.14.1
+julia 0.7
+StatsBase 0.24.0
+Distributions 0.16.1
 SpecialFunctions 0.3.8

docs/make.jl

@@ -24,7 +24,7 @@ makedocs(
 
 deploydocs(
     repo = "github.com/juliangehring/MultipleTesting.jl.git",
-    julia = "0.6",
+    julia = "0.7",
     target = "build",
     deps   = nothing,
     make   = nothing

src/MultipleTesting.jl

@@ -17,6 +17,8 @@ import Distributions: estimate
 
 import SpecialFunctions: digamma
 
+import Random: shuffle!
+
 export
     PValues,
     adjust,

src/combinations.jl

@@ -166,7 +166,7 @@ function combine(pValues::PValues{T}, method::StoufferCombination)::T where T<:A
     if minimum(pValues) == 0 || maximum(pValues) == 1
         return NaN
     end
-    z = cquantile.(Normal(), pValues)
+    z = cquantile.(Ref(Normal()), pValues)
     z = sum(z) ./ sqrt(n)
     p = ccdf(Normal(), z)
     return p
@@ -180,7 +180,7 @@ function combine(pValues::PValues{T}, weights::AbstractVector{R}, method::Stouff
     if minimum(pValues) == 0 || maximum(pValues) == 1
         return NaN
     end
-    z = cquantile.(Normal(), pValues) .* weights
+    z = cquantile.(Ref(Normal()), pValues) .* weights
     z = sum(z) ./ sqrt(sum(abs2, weights))
     p = ccdf(Normal(), z)
     return p
@@ -307,7 +307,7 @@ function combine(pValues::PValues{T}, method::WilkinsonCombination)::T where T<:
     if rank < 1 || rank > n
         throw(ArgumentError("Rank must be in 1,..,$(n)"))
     end
-    p_rank = select(pValues, rank) # faster than `sort(pValues)[rank]`
+    p_rank = partialsort(pValues, rank)
     p = cdf(Beta(rank, n-rank+1), p_rank)
     return p
 end

src/higher-criticism.jl

@@ -36,7 +36,7 @@ function estimate(pValues::PValues{T}, method::HigherCriticismScores) where T<:A
     denom = F .* (one(T) .- F) ./ n
     # avoid denominator of 0 for last value
     idx0 = denom .== 0
-    denom[idx0] = minimum(denom[.!idx0]) .+ eps()  # conservative
+    denom[idx0] .= minimum(denom[.!idx0]) .+ eps()  # conservative
     hcs = abs.(F .- pValues) ./ sqrt.(denom)
     return hcs
 end
@@ -68,6 +68,6 @@ struct HigherCriticismThreshold
 end
 
 function estimate(pValues::PValues{T}, method::HigherCriticismThreshold) where T<:AbstractFloat
-    idx_hcv = indmax(estimate(pValues, HigherCriticismScores()))
+    idx_hcv = argmax(estimate(pValues, HigherCriticismScores()))
     return pValues[idx_hcv]
 end

src/model.jl

@@ -43,8 +43,7 @@ partitioning the empirical distribution of p-values. Bioinformatics 19,
 
 """
 function BetaUniformMixtureModel(π0::AbstractFloat, α::AbstractFloat = 0.5, β::AbstractFloat = 1.0)
-    if !isin(π0, 0., 1.)
-        throw(DomainError())
-    end
-    MixtureModel([Beta(α, β), Uniform()], [1-π0, π0])
+    isin(π0, 0, 1) || throw(DomainError("π0 must be in [0, 1]"))
+    bum = MixtureModel([Beta(α, β), Uniform()], [1-π0, π0])
+    return bum
 end

src/pi0-estimators.jl

@@ -84,7 +84,10 @@ Series B (Statistical Methodology) 66, 187–205.
 struct Storey <: Pi0Estimator
     λ::Float64
 
-    Storey(λ) = isin(λ, 0, 1) ? new(λ) : throw(DomainError())
+    function Storey(λ)
+        isin(λ, 0, 1) || throw(DomainError("λ must be in [0, 1]"))
+        return new(λ)
+    end
 end
 
 Storey() = Storey(0.1)
@@ -131,8 +134,11 @@ struct StoreyBootstrap <: Pi0Estimator
     λseq::Vector{Float64}
     q   ::Float64
 
-    StoreyBootstrap(λseq, q) =
-        isin(λseq, 0, 1) && isin(q, 0, 1) ? new(λseq, q) : throw(DomainError())
+    function StoreyBootstrap(λseq, q)
+        isin(λseq, 0, 1) || throw(DomainError("λseq must be in [0, 1]"))
+        isin(q, 0, 1) || throw(DomainError("q must be in [0, 1]"))
+        return new(λseq, q)
+    end
 end
 
 StoreyBootstrap() = StoreyBootstrap(0.05:0.05:0.95, 0.1)
@@ -145,7 +151,7 @@ function estimate_pi0(pValues::PValues{T}, pi0estimator::StoreyBootstrap) where
     pi0 = w ./ n ./ (1 .- lambdas)
     min_pi0 = quantile(pi0, q)
     mse = (w ./ (n.^2 .* (1 .- lambdas).^2 )) .* (1 .- w/n) + (pi0 .- min_pi0).^2
-    pi0 = clamp(pi0[indmin(mse)], 0, 1)
+    pi0 = clamp(pi0[argmin(mse)], 0, 1)
     return pi0
 end
 
@@ -205,7 +211,7 @@ function lsl_pi0_vec(pValues::AbstractVector{T}) where T<:AbstractFloat
     pValues = sort_if_needed(pValues)
     s = (1 .- pValues) ./ (n:-1:1)
     d = diff(s) .< 0
-    idx = findfirst(d) + 1
+    idx = something(findfirst(d), 0) + 1
     pi0 = min( 1/s[idx] + 1, n ) / n
     return pi0
 end
@@ -229,7 +235,10 @@ julia> estimate_pi0(pvals, Oracle(0.5)) # a bit boring...
 struct Oracle <: Pi0Estimator
     π0::Float64
 
-    Oracle(π0) = isin(π0, 0, 1) ? new(π0) : throw(DomainError())
+    function Oracle(π0)
+        isin(π0, 0, 1) || throw(DomainError("π0 must be in [0, 1]"))
+        return new(π0)
+    end
 end
 
 Oracle() = Oracle(1.0)
@@ -269,7 +278,10 @@ struct TwoStep <: Pi0Estimator
     α::Float64
     adjustment::PValueAdjustment
 
-    TwoStep(α, method) = isin(α, 0, 1) ? new(α, method) : throw(DomainError())
+    function TwoStep(α, method)
+        isin(α, 0, 1) || throw(DomainError("α must be in [0, 1]"))
+        return new(α, method)
+    end
 end
 
 TwoStep() = TwoStep(0.05)
@@ -313,8 +325,10 @@ Statistical Society: Series B (Statistical Methodology) 74, 163–182.
 struct RightBoundary <: Pi0Estimator
     λseq::Vector{Float64}
 
-    RightBoundary(λseq) =
-        isin(λseq, 0, 1) ? new(λseq) : throw(DomainError())
+    function RightBoundary(λseq)
+        isin(λseq, 0, 1) || throw(DomainError("λseq must be in [0, 1]"))
+        return new(λseq)
+    end
 end
 
 # λseq used in Liang, Nettleton 2012
@@ -330,7 +344,7 @@ function estimate_pi0(pValues::PValues{T}, pi0estimator::RightBoundary) where T<
     pi0_estimates = reverse(cumsum(reverse(h.weights)))./(1 .- λseq)./n
     pi0_decrease = diff(pi0_estimates) .>= 0
     pi0_decrease[end] = true
-    pi0 = pi0_estimates[findfirst(pi0_decrease, true) + 1]
+    pi0 = pi0_estimates[something(findfirst(pi0_decrease), 0) + 1]
     pi0 = clamp(pi0, 0, 1)
     return pi0
 end
@@ -366,11 +380,11 @@ struct CensoredBUM <: Pi0Estimator
     maxiter::Int64
 
     function CensoredBUM(γ0, λ, xtol, maxiter)
-        if isin(γ0, 0, 1) && isin(λ, 0, 1) && isin(xtol, 0, 1) && maxiter > 0
-            new(γ0, λ, xtol, maxiter)
-        else
-            throw(DomainError())
-        end
+        isin(γ0, 0, 1) || throw(DomainError("γ0 must be in [0, 1]"))
+        isin(λ, 0, 1) || throw(DomainError("λ must be in [0, 1]"))
+        isin(xtol, 0, 1) || throw(DomainError("xtol must be in [0, 1]"))
+        maxiter > 0 || throw(DomainError("maxiter must be a positive number"))
+        return new(γ0, λ, xtol, maxiter)
     end
 end
 
@@ -451,7 +465,7 @@ function cbum_pi0_naive(pValues::AbstractVector{T},
         α = -sum(z[idx_right])
         α = α / ( ll * sum(z[idx_left]) + sum(z[idx_right] .* lpr) )
         xl = (1-γ) * (λ^α)
-        z[idx_left] = xl ./ (γ*λ + xl)
+        z[idx_left] .= xl / (γ*λ + xl)
         xr = (1-γ) * α * pValues[idx_right].^(α-1)
         z[idx_right] = xr ./ (γ .+ xr)
         pi0_new = γ + (1-γ)*α
@@ -495,11 +509,10 @@ struct BUM <: Pi0Estimator
     maxiter::Int64
 
     function BUM(γ0, xtol, maxiter)
-        if isin(γ0, 0, 1) && isin(xtol, 0, 1)
-            new(γ0, xtol, maxiter)
-        else
-            throw(DomainError())
-        end
+        isin(γ0, 0, 1) || throw(DomainError("γ0 must be in [0, 1]"))
+        isin(xtol, 0, 1) || throw(DomainError("xtol must be in [0, 1]"))
+        maxiter > 0 || throw(DomainError("maxiter must be a positive number"))
+        return new(γ0, xtol, maxiter)
     end
 end
 
@@ -620,11 +633,10 @@ struct ConvexDecreasing <: Pi0Estimator
     maxiter::Int64
 
     function ConvexDecreasing(gridsize, xtol, maxiter)
-        if gridsize > 0 && isin(xtol, 0, 1) && maxiter > 0
-            new(gridsize, xtol, maxiter)
-        else
-            throw(DomainError())
-        end
+        gridsize > 0 || throw(DomainError("gridsize must be a positive number"))
+        isin(xtol, 0, 1) || throw(DomainError("xtol must be in [0, 1]"))
+        maxiter > 0 || throw(DomainError("maxiter must be a positive number"))
+        return new(gridsize, xtol, maxiter)
     end
 end
 
@@ -713,11 +725,11 @@ function convex_decreasing(pValues::AbstractVector{T},
 end
 
 function find_theta(t::Vector{Float64}, p::Vector{Float64})
-    return indmax( [theta.^-2 * sum(theta .- p[p .< theta]) for theta in t] )
+    return argmax( [theta.^-2 * sum(theta .- p[p .< theta]) for theta in t] )
 end
 
 function find_theta(t::Vector{Float64}, p::Vector{Float64}, f_p::Vector{Float64})
-    return indmax( [theta.^-2 * sum( (theta .- p) .* (p .< theta) ./ f_p ) for theta in t] )
+    return argmax( [theta.^-2 * sum( (theta .- p) .* (p .< theta) ./ f_p ) for theta in t] )
 end
 
 function decide(f_p::Vector{Float64}, f_theta_p::Vector{Float64}, ε::Float64)

src/pval-adjustment.jl

@@ -726,4 +726,4 @@ function check_number_tests(k::Integer, n::Integer)
 end
 
 
-harmonic_number(n::Integer) =  digamma(n+1) + γ
+harmonic_number(n::Integer) =  digamma(n+1) + MathConstants.γ