# How to get the cartesian indexes of the k largest elements of a matrix (2D array) efficiently?

**URL:** <https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419>\
**Category:** New to Julia\
**Tags:** indexing, array\
**Created:** [October 8, 2021, 7:41am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419 "2021-10-08T07:41:56Z")\
**Posts on this page:** 13\
**Page:** 1

<div class="post-metadata">

**Author:** ![erwanlecarpentier](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/erwanlecarpentier/32/23448_2.png) [@erwanlecarpentier](https://discourse.julialang.org/u/erwanlecarpentier)\
**Post date:** [October 8, 2021, 7:41am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/1 "2021-10-08T07:41:56Z")

</div>

Hi! I would like to find a way to get the **cartesian** indexes of the `k` largest elements in a matrix (2D array). Do you have suggestions of a computationally efficient method to achieve that?

Example: for `k=3` with the matrix `A` below, the cartesian indexes would be `[(1, 3), (2, 2), (2, 3)]`.

```julia
julia> A = [1 1 3; 1 3 3]
2×3 Matrix{Int64}:
 1 1 3
 1 3 3

```

Thanks!

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [October 8, 2021, 9:29am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/2 "2021-10-08T09:29:30Z")

</div>

What Is the result expected for `k=2`?

---

<div class="post-metadata">

**Author:** ![erwanlecarpentier](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/erwanlecarpentier/32/23448_2.png) [@erwanlecarpentier](https://discourse.julialang.org/u/erwanlecarpentier)\
**Post date:** [October 8, 2021, 9:34am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/3 "2021-10-08T09:34:05Z")

</div>

I guess either ties are broken arbitrarily, or the first element is returned (I am interested in both solutions). By first element, I mean the first reached when going through the **linear indexing** of the matrix `A` (from 1 to 6 in this example).

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [October 8, 2021, 9:44am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/4 "2021-10-08T09:44:58Z")

</div>

Maybe something like

```julia
c = Vector{CartesianIndex{2}}(undef,k)
for i in lastindex(M):lastindex(M)-k
    c[i] = CartesianIndex...
end

```

(Sorry, on the cell phone…)

---

<div class="post-metadata">

**Author:** ![erwanlecarpentier](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/erwanlecarpentier/32/23448_2.png) [@erwanlecarpentier](https://discourse.julialang.org/u/erwanlecarpentier)\
**Post date:** [October 8, 2021, 10:19am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/5 "2021-10-08T10:19:04Z")

</div>

I am unsure I understood, but this made me think of 1) flatten the array, 2) sortperm to get the largest elements and 3) retrieve the cartesian indexes from the `k` last elements? This would give something like:

```julia
function klargest_indexes(m, k)
    ci = CartesianIndices(size(m))
    p = sortperm(vec(m))[end-k+1:end]
    ci[p]
end
k = 2
m = [1 2 3; 4 5 5; 5 4 3]
out = klargest_indexes(m, k)

```

In such a case, ties yield the last elements in the linear indexing. Is this what you hinted? Is it computationally efficient? I get a `@btime` of `5.374 μs (6 allocations: 3.55 KiB)` for `m = rand(20, 20)`, `k = 2`.

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [October 8, 2021, 10:54am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/6 "2021-10-08T10:54:52Z")

</div>

> [@erwanlecarpentier](#):
>
> `A = [1 1 3; 1 3 3]`

I meant probably something like this:

```julia
julia> ci = CartesianIndices(size(A))
       for i in lastindex(A)-2:lastindex(A) 
           @show ci[i]
       end
ci[i] = CartesianIndex(2, 2)
ci[i] = CartesianIndex(1, 3)
ci[i] = CartesianIndex(2, 3)

```

(you can by default run over a matrix as if it was a linearly indexed)

---

<div class="post-metadata">

**Author:** ![oxinabox](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oxinabox/32/206603_2.png) [@oxinabox](https://discourse.julialang.org/u/oxinabox)\
**Post date:** [October 8, 2021, 11:05am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/7 "2021-10-08T11:05:23Z")

</div>

> [@erwanlecarpentier](#):
>
> `sortperm(vec(m))[end-k+1:end]`

use instead

```julia
sortperm(vec(m), 1:k; rev=true)

```

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [October 8, 2021, 11:29am UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/8 "2021-10-08T11:29:12Z")

</div>

> [@oxinabox](#):
>
> use instead
> 
> ```julia
> sortperm(vec(m), 1:k; rev=true)
> 
> ```

I think you mean `partialsortperm(vec(m), 1:k; rev=true)`

---

<div class="post-metadata">

**Author:** ![erwanlecarpentier](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/erwanlecarpentier/32/23448_2.png) [@erwanlecarpentier](https://discourse.julialang.org/u/erwanlecarpentier)\
**Post date:** [October 8, 2021, 1:16pm UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/9 "2021-10-08T13:16:42Z")

</div>

> [@stevengj](#):
>
> partialsortperm(vec(m), 1:k; rev=true)

Thanks! This is indeed more efficient!

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [October 8, 2021, 2:06pm UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/10 "2021-10-08T14:06:28Z")

</div>

Am I missing something, or this is much more efficient?

```julia
julia> A = [1 1 3; 1 3 3]
2×3 Matrix{Int64}:
 1 1 3
 1 3 3

julia> function getci(M,k)
           ci = CartesianIndices(size(M))
           cis = Vector{CartesianIndex{2}}(undef,k)
           j = 0
           for i in lastindex(M)-k+1:lastindex(M)
               j += 1
               cis[j] = ci[i]
           end
           return cis
       end
getci (generic function with 1 method)

julia> @btime getci($A,3)
  32.426 ns (1 allocation: 128 bytes)
3-element Vector{CartesianIndex{2}}:
 CartesianIndex(2, 2)
 CartesianIndex(1, 3)
 CartesianIndex(2, 3)

```

vs

```julia
julia> function klargest_indexes(m, k)
           ci = CartesianIndices(size(m))
           p = partialsortperm(vec(m), 1:k; rev=true)
           ci[p]
       end
klargest_indexes (generic function with 1 method)

julia> @btime klargest_indexes($A,3)
  329.935 ns (6 allocations: 384 bytes)
3-element Vector{CartesianIndex{2}}:
 CartesianIndex(2, 2)
 CartesianIndex(1, 3)
 CartesianIndex(2, 3)

```

---

<div class="post-metadata">

**Author:** ![erwanlecarpentier](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/erwanlecarpentier/32/23448_2.png) [@erwanlecarpentier](https://discourse.julialang.org/u/erwanlecarpentier)\
**Post date:** [October 8, 2021, 2:37pm UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/11 "2021-10-08T14:37:13Z")

</div>

I think the sorting is lacking in the proposed `getci` function. Correct me if I am wrong but I believe it to yield the cartesian indexes of the `k` last elements rather than the `k` largest elements.

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [October 8, 2021, 2:38pm UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/12 "2021-10-08T14:38:58Z")

</div>

Ah, I see. In the example you posted those are the same, and I misunderstood your question. Thanks.

(I think you can do that without any sorting, just with an auxiliary array of size k, and running once over the elements of the matrix storing the indexes of the k larger)

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [October 8, 2021, 3:03pm UTC](https://discourse.julialang.org/t/how-to-get-the-cartesian-indexes-of-the-k-largest-elements-of-a-matrix-2d-array-efficiently/69419/13 "2021-10-08T15:03:20Z")

</div>

There you go:

```julia
julia> function getci!(M,k,cmins,vmins) # inplace
           ci = CartesianIndices(size(M))
           for i in 1:k
               cmins[i] = ci[i]
               vmins[i] = M[i]
           end
           imin = findmin(vmins)[2]
           for i in firstindex(M)+k:lastindex(M)
               if M[i] > vmins[imin]
                   cmins[imin] = ci[i]
                   vmins[imin] = M[i]
                   imin = findmin(vmins)[2]
               end
           end
           return cmins, vmins
       end

       function getci(M,k) # allocating
           cmins = Vector{CartesianIndex{2}}(undef,k)
           vmins = Vector{eltype(M)}(undef,k)
           return getci!(M,k,cmins,vmins)
       end
getci (generic function with 1 method)

julia> A = [3 1 1
            1 3 3]
2×3 Matrix{Int64}:
 3 1 1
 1 3 3

julia> @btime getci($A,3)
  78.857 ns (2 allocations: 240 bytes)
(CartesianIndex{2}[CartesianIndex(1, 1), CartesianIndex(2, 2), CartesianIndex(2, 3)], [3, 3, 3])

julia> cmins=Vector{CartesianIndex{2}}(undef,3); vmins=Vector{Int}(undef,3);

julia> @btime getci!($A,3,$cmins,$vmins)
  35.864 ns (0 allocations: 0 bytes)
(CartesianIndex{2}[CartesianIndex(1, 1), CartesianIndex(2, 2), CartesianIndex(2, 3)], [3, 3, 3])

```

(since here there is an auxiliary array with the associated values, why not return it)

(for `rand(20,20)` and `k=2`, an example you mentioned this is 10 times faster than the alternative:

```julia
julia> A = rand(20,20); k = 2;

julia> @btime klargest_indexes($A,2)
  5.540 μs (8 allocations: 3.56 KiB)
2-element Vector{CartesianIndex{2}}:
 CartesianIndex(14, 8)
 CartesianIndex(10, 10)

julia> @btime getci($A,2)
  549.492 ns (2 allocations: 208 bytes)
(CartesianIndex{2}[CartesianIndex(10, 10), CartesianIndex(14, 8)], [0.9941985771967297, 0.9992607929792117])

julia> cmins=Vector{CartesianIndex{2}}(undef,k); vmins=Vector{eltype(A)}(undef,k);

julia> @btime getci!($A,2,$cmins,$vmins);
  509.870 ns (0 allocations: 0 bytes)

```