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.

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

remove psignrank from Rmath #330

Open
wants to merge 2 commits into
base: master
Choose a base branch
from

Conversation

ArnoStrouwen
Copy link

@ArnoStrouwen ArnoStrouwen commented Feb 9, 2025

This PR attempts to reimplement psignrank in pure Julia, to get rid of any GNU license issues with the R version @andreasnoack .
The Rmath version returns the exact correct decimal number, while my version has some slight floating point error.
I'm not sure why this is because, the Rmath version also containts floating point operations:
https://github.com/JuliaStats/Rmath-julia/blob/6f2d37ff112914d65559bc3e0035b325c11cf361/src/signrank.c#L161

julia> @btime HypothesisTests.psignrank.(0:sum(0:n), n, true)
  1.030 μs (38 allocations: 5.41 KiB)
11-element Vector{Float64}:
 0.0625
 0.12500000000000003
 0.18750000000000003
 0.3125
 0.4375
 0.5625000000000001
 0.6875000000000001
 0.8125000000000001
 0.8749999999999999
 0.9375000000000001
 1.0
 
 julia> @btime Rmath.psignrank.(0:sum(0:n), n, true)
  1.950 μs (5 allocations: 256 bytes)
11-element Vector{Float64}:
 0.0625
 0.125
 0.1875
 0.3125
 0.4375
 0.5625
 0.6875
 0.8125
 0.875
 0.9375
 1.0


julia> @btime HypothesisTests.psignrank.(0:sum(0:n), n, false)
  1.130 μs (38 allocations: 5.41 KiB)
11-element Vector{Float64}:
 1.0
 0.9375000000000001
 0.8749999999999999
 0.8125000000000001
 0.6875000000000001
 0.5625000000000001
 0.4375
 0.3125
 0.18750000000000003
 0.12500000000000003
 0.0625


julia> @btime Rmath.psignrank.((0:sum(0:n)).-1, n, false)
  1.960 μs (6 allocations: 288 bytes)
11-element Vector{Float64}:
 1.0
 0.9375
 0.875
 0.8125
 0.6875
 0.5625
 0.4375
 0.3125
 0.1875
 0.125
 0.0625

@codecov-commenter
Copy link

codecov-commenter commented Feb 9, 2025

Codecov Report

Attention: Patch coverage is 80.00000% with 4 lines in your changes missing coverage. Please review.

Project coverage is 93.49%. Comparing base (1cbeccd) to head (739d3b5).

Files with missing lines Patch % Lines
src/wilcoxon.jl 80.00% 4 Missing ⚠️
Additional details and impacted files
@@            Coverage Diff             @@
##           master     #330      +/-   ##
==========================================
- Coverage   94.04%   93.49%   -0.55%     
==========================================
  Files          29       29              
  Lines        1814     1830      +16     
==========================================
+ Hits         1706     1711       +5     
- Misses        108      119      +11     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

@ArnoStrouwen
Copy link
Author

ArnoStrouwen commented Feb 9, 2025

Julia 1.3 complains about the use of begin, working around would be not too difficult, or we can up the CI bound to something more modern?

Copy link
Member

@devmotion devmotion left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I wonder if the function should be added to StatsFuns, to make it easier and to reduce dependencies in case you want to use it in some other package and since it's currently the main package with Julia implementations of Rmath functions.

src/wilcoxon.jl Outdated Show resolved Hide resolved
src/wilcoxon.jl Outdated Show resolved Hide resolved
src/wilcoxon.jl Outdated Show resolved Hide resolved
src/wilcoxon.jl Outdated Show resolved Hide resolved
src/wilcoxon.jl Outdated
int_type = typeof(n)
W = int_type(W)
two = one(int_type) + one(int_type)
max_W = n * (n + one(int_type)) ÷ two + one(int_type)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is it actually possible to construct a Vector with non-Int length? I assume this could just always be an Int as it doesn't have any ilpact on the calculated values.

src/wilcoxon.jl Outdated Show resolved Hide resolved
src/wilcoxon.jl Outdated Show resolved Hide resolved
@ArnoStrouwen
Copy link
Author

I wonder if the function should be added to StatsFuns, to make it easier and to reduce dependencies in case you want to use it in some other package and since it's currently the main package with Julia implementations of Rmath functions.

It should really go to Distributions, since this is calculating the null distribution, but I'm too lazy to implement all the other functions. I would leave it here for now and create an issue for upstreaming proper Wilcoxon related distributions.

@ArnoStrouwen
Copy link
Author

Not bad:

julia> n = 4
       @btime HypothesisTests.psignrank.(0:sum(0:n), n, true)
  401.493 ns (25 allocations: 1.08 KiB)
11-element Vector{Float64}:
 0.0625
 0.125
 0.1875
 0.3125
 0.4375
 0.5625
 0.6875
 0.8125
 0.875
 0.9375
 1.0

julia> @btime HypothesisTests.psignrank.(0:sum(0:n), n, false)
  422.111 ns (25 allocations: 1.08 KiB)
11-element Vector{Float64}:
 1.0
 0.9375
 0.875
 0.8125
 0.6875
 0.5625
 0.4375
 0.3125
 0.1875
 0.125
 0.0625

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants