From 536af3ef45a16750b8ad07514f0d8a9fdb8e5c86 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Kristian=20Agas=C3=B8ster=20Haaga?= Date: Mon, 5 Aug 2024 15:38:37 +0200 Subject: [PATCH] Fix implementation for mu < 0 --- src/independence_tests/secmi/secmi_test.jl | 22 ++++++++++++++++------ 1 file changed, 16 insertions(+), 6 deletions(-) diff --git a/src/independence_tests/secmi/secmi_test.jl b/src/independence_tests/secmi/secmi_test.jl index 6da990f2..2e4af9c2 100644 --- a/src/independence_tests/secmi/secmi_test.jl +++ b/src/independence_tests/secmi/secmi_test.jl @@ -89,16 +89,26 @@ function independence(test::SECMITest, x, y, z) σ̂ = 1/(nshuffles - 1) * sum((sₖ - μ̂)^2 for sₖ in secmiₖ) emp_cdf = ecdf(secmiₖ) F𝒩 = Normal(μ̂, σ̂) - # Degrees of freedom for Chi squared distribution estimated as the mean of the `secmiₖ` - # (page 18 in Kubkowski et al.). - F𝒳² = Chisq(μ̂) - D𝒩, D𝒳² = sup_values(emp_cdf, F𝒩, F𝒳², secmiₖ) - if D𝒩 < D𝒳² || μ̂ ≤ 1.0 + + if μ̂ ≤ 0.0 p = 1 - cdf(F𝒩, secmi₀) + return SECMITestResult(3, secmi₀, secmiₖ, p, μ̂, σ̂, emp_cdf, nothing, nothing) else - p = 1 - cdf(F𝒳², secmi₀) + # Degrees of freedom for Chi squared distribution estimated as the mean of the `secmiₖ` + # (page 18 in Kubkowski et al.). The `Chisq` distribution is only defined for μ̂ > 0, + # so we put μ̂ <= 0.0 in a separate criterion first to avoid errors. + F𝒳² = Chisq(μ̂) + D𝒩, D𝒳² = sup_values(emp_cdf, F𝒩, F𝒳², secmiₖ) + if D𝒩 < D𝒳² + p = 1 - cdf(F𝒩, secmi₀) + else + p = 1 - cdf(F𝒳², secmi₀) + end + return SECMITestResult(3, secmi₀, secmiₖ, p, μ̂, σ̂, emp_cdf, D𝒩, D𝒳²) + end + return SECMITestResult(3, secmi₀, secmiₖ, p, μ̂, σ̂, emp_cdf, D𝒩, D𝒳²) end