# How to speed up rowsum function?

**URL:** <https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664>\
**Category:** Performance\
**Tags:** question\
**Created:** [October 13, 2022, 1:19pm UTC](https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664 "2022-10-13T13:19:12Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![Strange\_Xue](https://avatars.discourse-cdn.com/v4/letter/s/e9a140/32.png) [@Strange\_Xue](https://discourse.julialang.org/u/Strange_Xue)\
**Post date:** [October 13, 2022, 1:19pm UTC](https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664/1 "2022-10-13T13:19:12Z")

</div>

There is a `rowsum` function in R, it’s very helpful and fast when constructing some likelihood function, `rowsum` can apply a function to a group subsetted from a matrix then concatenate these resulted vectors to a new matrix. So I write a similar funtion myself in Julia, however, the speed is no ideal, julia version costs 500ms, R version costs 5ms. Is there any way to speed up julia version code ?

```julia
# Julia version
using Random
using BenchmarkTools
seed = Random.seed!(2022 - 20 - 13)

function rowsum(m::AbstractMatrix{<:Number}, g::Vector{Int})
    rst = Matrix{eltype(m)}(undef, maximum(g), size(m, 2))
    @inbounds @simd for group in 1:maximum(g)
        rst[group, :] = sum(view(m, g .== group, :); dims = 1)
    end 
    return rst
end

N = 100000
a = rand(seed, N, 4);
g = repeat(1:Int(N/5), inner = 5);

rowsum(a, g);

@btime rowsum($a, $g);
@code_warntype rowsum(a, g) # make sure type stability

```

```nohighlight
# R version
library(microbenchmark)
set.seed(2022 - 20 - 13)

N <- 100000
a <- matrix(rnorm(N * 4), N, 4)

g <- rep(1:(N/5), each = 5)

b <- rowsum(a, g)

t = microbenchmark(
    "rowsum" = {
        b = rowsum(a, g)
    }, 
    times = 1000,
    unit = "ms"
)

t

```

---

<div class="post-metadata">

**Author:** ![oheil](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oheil/32/220745_2.png) [@oheil](https://discourse.julialang.org/u/oheil)\
**Post date:** [October 13, 2022, 1:58pm UTC](https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664/2 "2022-10-13T13:58:37Z")

</div>

Your Julia implementation is not efficient.  
Your grouping in

```julia
view(m, g .== group, :)

```

has to check the complete vector `g` again and again for every value in `group`.  
Per definition of the R function the same can be done on a single run through `g`, e.g. by storing the row index for each factor found in `g` in a dictionary.

In other words, your implementation compares each `group` value `maximum(g) * 100000` times, means

```julia
julia> length(1:maximum(g))*length(g)
2000000000

```

This times comparisons vs.

```julia
julia> length(g)
100000

```

this times running a single run through g.

I am not sure if I am clear enough. In the mean time I will try some more efficient implementation.

---

<div class="post-metadata">

**Author:** ![mcabbott](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcabbott/32/6603_2.png) [@mcabbott](https://discourse.julialang.org/u/mcabbott)\
**Post date:** [October 13, 2022, 2:17pm UTC](https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664/3 "2022-10-13T14:17:26Z")

</div>

Here’s one way to improve this, but there might be better ways:

```julia
julia> function rowsum2(mat::AbstractMatrix{<:Number}, groups::AbstractVector{<:Integer})
           size(mat, 1) == length(groups) || error("length of group vector must match matrix")
           rst = fill!(similar(mat, maximum(groups), size(mat, 2)), 0)
           for (i, g) in pairs(groups)
             @views rst[g, :] .+= mat[i, :]
           end
           return rst
       end;

julia> r = @btime rowsum($a, $g);
  min 1.125 s, mean 1.127 s (160002 allocations, 331.12 MiB)

julia> r ≈ @btime rowsum2($a, $g)
  min 984.500 μs, mean 1.061 ms (2 allocations, 625.05 KiB)
true

```

---

<div class="post-metadata">

**Author:** ![oheil](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oheil/32/220745_2.png) [@oheil](https://discourse.julialang.org/u/oheil)\
**Post date:** [October 13, 2022, 2:21pm UTC](https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664/4 "2022-10-13T14:21:16Z")

</div>

Here is my implementation, which isn’t faster but perhaps for OP better to read and which matches better to my wall of text from above:

```julia
function rowsum3(m::AbstractMatrix{<:Number}, g::Vector{Int})
	rst = Matrix{eltype(m)}(undef, maximum(g), size(m, 2))
	idx = [Int[] for i in 1:maximum(g) ]
	for i in 1:length(g)
       push!(idx[g[i]],i)
    end
    @inbounds for (i,rows) in enumerate(idx)
        rst[i, :] = sum(view(m, rows, :); dims = 1)
    end 
	return rst
end

```

```julia
julia> @btime rowsum($a, $g);
  965.716 ms (160002 allocations: 331.12 MiB)

julia> @btime rowsum3($a, $g);
  7.278 ms (120004 allocations: 8.09 MiB)

```

(EDIT: gave it the name rowsum3)

---

<div class="post-metadata">

**Author:** ![Strange\_Xue](https://avatars.discourse-cdn.com/v4/letter/s/e9a140/32.png) [@Strange\_Xue](https://discourse.julialang.org/u/Strange_Xue)\
**Post date:** [October 13, 2022, 3:52pm UTC](https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664/5 "2022-10-13T15:52:10Z")

</div>

You’re right. It’s really a big problem

---

<div class="post-metadata">

**Author:** ![Strange\_Xue](https://avatars.discourse-cdn.com/v4/letter/s/e9a140/32.png) [@Strange\_Xue](https://discourse.julialang.org/u/Strange_Xue)\
**Post date:** [October 13, 2022, 3:52pm UTC](https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664/6 "2022-10-13T15:52:58Z")

</div>

It’s much faster than R version, great solution!

---

<div class="post-metadata">

**Author:** ![rocco\_sprmnt21](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rocco_sprmnt21/32/20127_2.png) [@rocco\_sprmnt21](https://discourse.julialang.org/u/rocco_sprmnt21)\
**Post date:** [October 13, 2022, 4:57pm UTC](https://discourse.julialang.org/t/how-to-speed-up-rowsum-function/88664/7 "2022-10-13T16:57:40Z")

</div>

A slight modification to the @mcabbott solution that makes it a little faster (running the matrices by column), but perhaps loses too much in generality.  
I do not know if the use cases where you have to apply it do not fall within the scheme that I adopted (groups of fixed length, in essence)

```julia
julia> function rowsum3(mat::AbstractMatrix{<:Number}, groups::AbstractVector{<:Integer})
                  size(mat, 1) == length(groups) || error("length of group vector must match matrix")
                  rst = fill!(similar(mat, maximum(groups), size(mat, 2)), 0)    
                   gr = maximum(groups)*size(mat, 2)
                   sgr = Int(length(groups)/maximum(groups))
                   for g in 1:gr
                       for i in (g-1)*sgr+1:g*sgr
                           @views rst[g] += mat[i]
                       end
                   end
                  return rst
              end;

julia> @btime rowsum3($a, $g);
  552.300 μs (2 allocations: 625.05 KiB)

julia> rowsum2(a, g) == rowsum3(a, g) 
true

```