I have a vector A
. I need a vector RA
which ranks each element of A
. Is there a better way to do this than what I have below?
A = [4;2;5;1;6;7]
RA = sum((A .== sort(A)') .* cumsum(ones(length(A)))', dims=2)
display([A RA])
I have a vector A
. I need a vector RA
which ranks each element of A
. Is there a better way to do this than what I have below?
A = [4;2;5;1;6;7]
RA = sum((A .== sort(A)') .* cumsum(ones(length(A)))', dims=2)
display([A RA])
julia> sortperm(A)
6-element Vector{Int64}:
4
2
1
3
5
6
julia> A[sortperm(A)]
6-element Vector{Int64}:
1
2
4
5
6
7
do you want some kind of sortperm?
I think they mean the inverse of sortperm
, the solution would be invperm(sortperm(A))
.
ehh, sortperm(; rev = true)
?
If I understand the problem correctly, the result should be vector with indices of the elements after sorting, for A
it would be [ 3 2 4 1 5 6]
. Command sortperm(A, rev = true)
does not return it, invperm(sortperm(A))
does.
Thanks invperm(sortperm(A))
solves it
just as an alternative to the first proposal
[findfirst(sa->sa==a, sort(A)) for a in A]
This function built on a simple loop, could be an extra argument in favor of the thesis (which is not mine: I prefer the invperm (...)
) solution) that sometimes a for loop may be preferable to the use of Base library functions
function rankfy1(A)
rank=Int[]
for a in A
r=1
for i in eachindex(A)
if a > A[i]
r+=1
end
end
push!(rank, r)
end
return rank
end
#or better
function rankfy2(A)
rank=similar(A)
for i in eachindex(A)
r=1
for ii in eachindex(A)
if A[i] > A[ii]
r+=1
end
end
rank[i]=r
end
return rank
end
If we are playing with the code this algorithm can be also expressed as [count(<=(a),A) for a in A]
; if we do not want Base
one can also write the same thing like this
map(A) do a
S = 0
for b in A
b <= a && (S += 1)
end
S
end
All are O(n^2), so not very good.
OK.
I had done some quick tests on the vector A proposed by the OP.
On larger vectors invperm () proves to be much more efficient.
This confirms me in the idea that, whenever possible, it is better to use functions written by “expert” developers rather than relying on your own loops.
This, for example, does not perform as invperm (sortperm ())
but is much better than just looping.
A=rand(1:10^8, 10^5)
last.(sort(tuple.(sort(tuple.(A,1:length(A)), by=first), 1:length(A)), by=last∘first))
PS
what is the non-O (n ^ 2) algorithm used by invperm (…)?
See Rankings and Rank Correlations · StatsBase.jl for general ranking: it lets you disambiguate ranks when ties are present.
It seems to do something like this, obviously in a much more general way.
function rank1(A)
sp=sortperm(A)
invsp=Array{Int}(undef, length(A))
foreach(i->invsp[sp[i]]=i, eachindex(A))
return invsp
end
In Julia you can check where the body of the function is with @which
and even display it directly by @edit
. If you use them you will get invperm
algorithm. This algorithm is essentially simple
function myInvperm(p)
res = Vector{Int}(undef,size(p))
res[p] .= 1:length(p)
end
This is O(n)
Thank you.
Especially for the last line
res [p] = 1: length (p)
This allows me to rewrite the function in the following way.
PS
I noticed that explicitly using the return invsp, @btime
calculates one less allocation.
julia> function rank2(A)
sp=sortperm(A)
invsp=Array{Int}(undef, length(A))
invsp[sp]=eachindex(A)
end
rank2 (generic function with 1 method)
julia> @btime rank2(A);
5.804 ms (6 allocations: 1.53 MiB)
julia> function rank2(A)
sp=sortperm(A)
invsp=Array{Int}(undef, length(A))
invsp[sp]=eachindex(A)
return invsp
end
rank2 (generic function with 1 method)
julia> @btime rank2(A);
5.834 ms (5 allocations: 1.53 MiB)