fix szcf-weiya/paperMonoDecomp#32

szcf-weiya · Oct 17, 2023 · 71653be · 71653be
1 parent e684fbc
commit 71653be
Show file tree

Hide file tree

Showing 4 changed files with 16 additions and 4 deletions.
diff --git a/src/mono_decomp.jl b/src/mono_decomp.jl
@@ -450,8 +450,12 @@ function cv_mono_decomp_ss(x::AbstractVector{T}, y::AbstractVector{T}; figname =
     s_smoothness = λ * (γup .+ γdown)' * Ω * (γup .+ γdown)
     # s_discrepancy = [max(eps(), min(s_residual, s_smoothness)) / 10^k_magnitude, 
     #                             max(s_residual, s_smoothness) * 10^k_magnitude]
-    s_discrepancy = [eps(), max(s_residual, s_smoothness) * 10^k_magnitude]
+    s_discrepancy = [eps(), max(s_residual, s_smoothness, eps()) * 10^k_magnitude]
+    if s0 == 0
+        μrange = [eps(), eps() * 10^k_magnitude]
+    else
         μrange = s_discrepancy / s0^2
+    end
     verbose && @info "μrange for compatible terms: $μrange"
     verbose && @info "s0 = $s0, s_residual = $s_residual, s_smoothness = $s_smoothness, s_discrepancy = $s_discrepancy"
     if μrange[2] < minmaxμ # cbrt of eps()

diff --git a/src/mono_test.jl b/src/mono_test.jl
@@ -311,7 +311,7 @@ function mono_test_bootstrap_cs(x::AbstractVector{T}, y::AbstractVector{T}; nrep
     end
     @debug ts
     @debug tobs
-    pval = sum(ts .> tobs) / nrep
+    pval = (sum(ts .> tobs) + sum(ts .== tobs) * 0.5) / nrep
     return pval, D
 end
 
@@ -537,7 +537,7 @@ function mono_test_bootstrap_ss(x::AbstractVector{T}, y::AbstractVector{T}; nrep
             ts[i] = Inf
         end
     end
-    pval = sum(ts .> tobs) / nrep
+    pval = (sum(ts .> tobs) + sum(ts .== tobs) * 0.5) / nrep
     return pval, D
 end
 

diff --git a/test/test_monodecomp.jl b/test/test_monodecomp.jl
@@ -428,7 +428,7 @@ end
 end
 
 @testset "GI solver for MDCS" begin
-    tol = 1e-5
+    tol = 1e-4
     x, y, _ = gen_data(100, 0.1, x->x^3)
     D1 = mono_decomp_cs(x, y, s = 1.0, J = 10, use_GI = true)
     D2 = mono_decomp_cs(x, y, s = 1.0, J = 10, use_GI = false)

diff --git a/test/test_monotest.jl b/test/test_monotest.jl
@@ -12,6 +12,14 @@ using MonotoneDecomposition
     @test MonotoneDecomposition.mono_test_bootstrap_ss(x, x.^2)[1] < 0.05
 end
 
+@testset "both components are zero" begin
+    x = randn(100)
+    y = zeros(100)
+    @test MonotoneDecomposition.mono_test_bootstrap_cs(x, x)[1] > 0.05
+    @test MonotoneDecomposition.mono_test_bootstrap_cs(x, x, h0_mono = true)[1] > 0.05
+    @test MonotoneDecomposition.mono_test_bootstrap_ss(x, x)[1] > 0.05
+end
+
 @testset "ghosal constant" begin
     @test MonotoneDecomposition.ghosal_λ() ≈ 9.974576271186443 atol=1e-5
     c = MonotoneDecomposition.ghosal_critical_value()