Ranking of elements of a vector

General Usage
,
1

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])
2 
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?

3

I think they mean the inverse of sortperm, the solution would be invperm(sortperm(A)).

1 Like
4

ehh, sortperm(; rev = true)?

5

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.

2 Likes
6

Thanks invperm(sortperm(A)) solves it

7

just as an alternative to the first proposal :grinning:

[findfirst(sa->sa==a, sort(A)) for a in A]
8

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
9

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.

2 Likes
10

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 (…)?

11

See Rankings and Rank Correlations · StatsBase.jl for general ranking: it lets you disambiguate ranks when ties are present.

1 Like
12

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
13

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)

1 Like
14

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)