# Multithreaded broadcast on a grouped DimensionalData array

**URL:** <https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607>\
**Category:** General Usage\
**Tags:** multithreading, broadcast, loopvectorization, dimensionaldata, rasters\
**Created:** [July 29, 2024, 6:11pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607 "2024-07-29T18:11:51Z")\
**Posts on this page:** 14\
**Page:** 2

<div class="post-metadata">

**Author:** ![Raf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raf/32/3383_2.png) [@Raf](https://discourse.julialang.org/u/Raf)\
**Post date:** [July 31, 2024, 6:54pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/21 "2024-07-31T18:54:20Z")

</div>

The groupby array is new and may not be fully tested for generic broadcast. A GitHub issue is the best place to resolve this.

(And sorry for the delay, I’m on holiday)

---

<div class="post-metadata">

**Author:** ![aasdelat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/aasdelat/32/45268_2.png) [@aasdelat](https://discourse.julialang.org/u/aasdelat)\
**Post date:** [July 31, 2024, 7:19pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/22 "2024-07-31T19:19:56Z")

</div>

Don’t worry, it is your duty to enjoy your holidays!

I will make the GitHub issue, but, please, don’t return to this while you are on holidays.

Have fun!

---

<div class="post-metadata">

**Author:** ![aasdelat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/aasdelat/32/45268_2.png) [@aasdelat](https://discourse.julialang.org/u/aasdelat)\
**Post date:** [August 3, 2024, 5:50pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/23 "2024-08-03T17:50:56Z")

</div>

> [@danielwe](#):
>
> ```julia
> using DimensionalData, Dates, Statistics, OhMyThreads
> 
> tempo = range(DateTime(2000); step=Hour(1), length=365*24*2)
> A = rand(X(1:0.01:2), Ti(tempo))
> groups = groupby(A, Ti => month)
> 
> tmean(g) = mean(g; dims=:Ti)
> tmap(tmean, groups)
> 
> ```

Hi, @danielwe , can you repeat your attempt, but after starting julia with the option --treads=# (# = number of desired threads)?

It seems that it does not work when multiple threads are enabled.

---

<div class="post-metadata">

**Author:** ![danielwe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielwe/32/35657_2.png) [@danielwe](https://discourse.julialang.org/u/danielwe)\
**Post date:** [August 3, 2024, 6:34pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/24 "2024-08-03T18:34:58Z")

</div>

You’re right, I get a `MethodError` for `resize!(::DimGroupByArray,...)`. Sorry about that, I thought I had used multiple threads when trying this out, but apparently not.

---

<div class="post-metadata">

**Author:** ![aasdelat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/aasdelat/32/45268_2.png) [@aasdelat](https://discourse.julialang.org/u/aasdelat)\
**Post date:** [August 3, 2024, 6:36pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/25 "2024-08-03T18:36:29Z")

</div>

So, could this be an issue for OhMyThreads instead of Rasters?

---

<div class="post-metadata">

**Author:** ![danielwe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielwe/32/35657_2.png) [@danielwe](https://discourse.julialang.org/u/danielwe)\
**Post date:** [August 4, 2024, 4:17am UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/26 "2024-08-04T04:17:12Z")

</div>

No, I wouldn’t say that. DimensionalData defines some special data types and implements dedicated methods for them. For example, `map(tmean, groups)` works because of a special method of `Base.map` implemented in DimensionalData, see here: [DimensionalData.jl/src/array/methods.jl at 2012cded9dab12cc1f622a082d4b22cdefeb1e7b · rafaqz/DimensionalData.jl · GitHub](https://github.com/rafaqz/DimensionalData.jl/blob/2012cded9dab12cc1f622a082d4b22cdefeb1e7b/src/array/methods.jl#L97-L101)

When Julia is running single-threaded, `tmap` from OhMyThreads falls back to `Base.map`, so everything works. When Julia is running with several threads, `tmap` tries to do its own multithreaded thing, but it doesn’t work because the types in DimensionalData need special treatment. Packages like OhMyThreads work through generic interfaces and can’t be held responsible for such special cases.

I’d like to emphasize that this is almost certainly a very simple problem for someone familiar with DimensionalData. Distributing the operations across threads is straightforward, the part that’s causing trouble is reassembling the results. Unfortunately, I’ve never used DimensionalData, I just responded because I’ve got some experience with LoopVectorization and multithreading. Hopefully, someone who knows DimensionalData well will chime in eventually.

---

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [August 4, 2024, 10:25am UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/27 "2024-08-04T10:25:24Z")

</div>

> [@danielwe](#):
>
> DimensionalData also seems to overload the broadcasting mechanism such that `mean.(datavar_grouped)` returns a `DimArray` instead of a plain array, as shown here: [Group By | DimensionalData.jl](https://rafaqz.github.io/DimensionalData.jl/v0.27.5/groupby#Grouping-and-reducing). I think this is why FastBroadcast falls back to regular broadcasting and refuses to apply any optimizations, including multithreading.

Yes, FastBroadcast only applies to broadcasts with `DefaultArrayStyle`. Otherwise, it falls back to whatever the custom style’s implementation is.

---

<div class="post-metadata">

**Author:** ![Raf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raf/32/3383_2.png) [@Raf](https://discourse.julialang.org/u/Raf)\
**Post date:** [August 4, 2024, 11:53am UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/28 "2024-08-04T11:53:15Z")

</div>

Hi all, sorry it took a while to get to this.

May I suggest just a threaded loop, before getting too complicated:

```julia
using DimensionalData, Dates, Statistics

tempo = range(DateTime(2000); step=Hour(1), length=365*24*2)
A = rand(X(1:0.01:2), Ti(tempo))
groups = groupby(A, Ti => month)

newdimarray = similar(groups, Float64)

```

Then:

```julia
julia> @btime for i in eachindex(newdimarray)
          newdimarray[i] = mean(groups[i])
       end
  5.837 ms (181 allocations: 12.22 KiB)

julia> @btime Threads.@threads :static for i in eachindex(newdimarray)
          newdimarray[i] = mean(groups[i])
       end
  2.963 ms (212 allocations: 16.02 KiB)

```

For the packages, it seems the main problems were `append!` and `resize!` not working.

There are two reasons this is difficult in DimensionalData.jl. The first is that irregular, unordered and categorical lookups cant be resized larger as we don’t know the new lookup values. Those are not actually problems here. The second is that a `DimArray` can hold a view and other lazy objects that also don’t define `append! or `resize!` - as with a DimGroupByArray here, it normally holds some other wrapper array that does the grouping. I guess we could just error in all those cases and work wherever possible.

(there is another more difficult reason. A `DimArray` has mostly immutable fields for performance and GPU compatability, so you can’t actually update an `AbstractRange` lookup in `resize!` or `append!` - and the lookup size always has to match the array size so we simply cant do it. Use-cases like this need to use `cat` rather than `append!` to just return a new object)

With broadcasts, `DimensionalStyle` is what lets a `DimArray` or a `Raster` keep its lookups and metadata through broadcasts, and is pretty fundamental to the package. If another package needs `DefaultArrayStyle` then either DimensionalData.jl or that package needs an extension with dispatch on `DimensionalStyle`.

---

<div class="post-metadata">

**Author:** ![aasdelat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/aasdelat/32/45268_2.png) [@aasdelat](https://discourse.julialang.org/u/aasdelat)\
**Post date:** [August 5, 2024, 6:46pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/29 "2024-08-05T18:46:46Z")

</div>

Thank you, @Raf .

It does not work. Each element of the outer `Ti` dimension is not a Float64, but another raster with two dimensions: `X` and `Ti`.

```julia
julia> groups
╭───────────────────────────────────────────────────╮
│ 12-element DimGroupByArray{DimArray{Float64,1},1} │
├───────────────────────────────────────────────────┴───────────── dims ┐
  ↓ Ti Sampled{Int64} [1, 2, …, 11, 12] ForwardOrdered Irregular Points
├───────────────────────────────────────────────────────────── metadata ┤
  Dict{Symbol, Any} with 1 entry:
  :groupby => :Ti=>month
├─────────────────────────────────────────────────────────── group dims ┤
  ↓ X, → Ti
└───────────────────────────────────────────────────────────────────────┘
  1 101×1488 DimArray
  2 101×1368 DimArray
  ⋮  
 11 101×1440 DimArray
 12 101×1464 DimArray

```

Your code gives:

```julia
julia> newdimarray
╭────────────────────────────────╮
│ 12-element DimArray{Float64,1} │
├────────────────────────────────┴──────────────────────────────── dims ┐
  ↓ Ti Sampled{Int64} [1, 2, …, 11, 12] ForwardOrdered Irregular Points
├───────────────────────────────────────────────────────────── metadata ┤
  Dict{Symbol, Any} with 1 entry:
  :groupby => :Ti=>month
└───────────────────────────────────────────────────────────────────────┘
  1 6.36837e-310
  2 6.36837e-310
  3 6.36837e-310
  ⋮  
 10 3.1e-322
 11 9.7e-322
 12 6.36838e-310

```

And, in the `mean()` call, we have to specify dims=:Ti.  
So, `newdimarray` should have the same structure as groups, but with the inner Rasters having dimensions `101×1` instead of `101×1488`.  
I have attempted to create such an structure, but could not figure out how.

---

<div class="post-metadata">

**Author:** ![Raf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raf/32/3383_2.png) [@Raf](https://discourse.julialang.org/u/Raf)\
**Post date:** [August 6, 2024, 9:44am UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/30 "2024-08-06T09:44:27Z")

</div>

It was an example… you can equally make a DimArray of Rasters. You can even just take the mean of the first one to get the type.

`out = similar(groups, typeof(mean(first(groups))); dims=Ti))`

Or just do `out = similar(groups, Raster)` to leave it loosely typed.

I get it that easier threading without defining the array would be nice here. But adding threading to DimensionalData.jl directly is pretty far out of scope, so to make it work someone just has to add an extension on a threading package. I would totally merge that, but I don’t have time or need to write it myself.

---

<div class="post-metadata">

**Author:** ![aasdelat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/aasdelat/32/45268_2.png) [@aasdelat](https://discourse.julialang.org/u/aasdelat)\
**Post date:** [August 6, 2024, 2:59pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/31 "2024-08-06T14:59:32Z")

</div>

Thank you, @Raf , it works!

In fact, I already tested this solution, but I got an error. I realize now that the error is with regards to showing the data, but not necessarily to the data itself (please, confirm). So, just adding a semicolon at the end of the line, avoids the error.

I dump here the error:

```julia
newdimarray = similar(groups, typeof(mean(first(groups), dims=Ti)) )
╭───────────────────────────────────────────────────╮
│ 12-element DimGroupByArray{DimArray{Float64,2},1} │
├───────────────────────────────────────────────────┴───────────── dims ┐
  ↓ Ti Sampled{Int64} [1, 2, …, 11, 12] ForwardOrdered Irregular Points
├───────────────────────────────────────────────────────────── metadata ┤
  Dict{Symbol, Any} with 1 entry:
  :groupby => :Ti=>month
Error showing value of type DimensionalData.DimGroupByArray{DimMatrix{Float64, Tuple{X{DimensionalData.Dimensions.Lookups.Sampled{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, DimensionalData.Dimensions.Lookups.ForwardOrdered, DimensionalData.Dimensions.Lookups.Regular{Float64}, DimensionalData.Dimensions.Lookups.Points, DimensionalData.Dimensions.Lookups.NoMetadata}}, Ti{DimensionalData.Dimensions.Lookups.Sampled{DateTime, Vector{DateTime}, DimensionalData.Dimensions.Lookups.ForwardOrdered, DimensionalData.Dimensions.Lookups.Irregular{Tuple{DateTime, DateTime}}, DimensionalData.Dimensions.Lookups.Points, DimensionalData.Dimensions.Lookups.NoMetadata}}}, Tuple{}, Matrix{Float64}, DimensionalData.NoName, DimensionalData.Dimensions.Lookups.NoMetadata}, 1, Tuple{Ti{DimensionalData.Dimensions.Lookups.Sampled{Int64, Vector{Int64}, DimensionalData.Dimensions.Lookups.ForwardOrdered, DimensionalData.Dimensions.Lookups.Irregular{Tuple{Nothing, Nothing}}, DimensionalData.Dimensions.Lookups.Points, DimensionalData.Dimensions.Lookups.NoMetadata}}}, Tuple{}, Vector{DimMatrix{Float64, Tuple{X{DimensionalData.Dimensions.Lookups.Sampled{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, DimensionalData.Dimensions.Lookups.ForwardOrdered, DimensionalData.Dimensions.Lookups.Regular{Float64}, DimensionalData.Dimensions.Lookups.Points, DimensionalData.Dimensions.Lookups.NoMetadata}}, Ti{DimensionalData.Dimensions.Lookups.Sampled{DateTime, Vector{DateTime}, DimensionalData.Dimensions.Lookups.ForwardOrdered, DimensionalData.Dimensions.Lookups.Irregular{Tuple{DateTime, DateTime}}, DimensionalData.Dimensions.Lookups.Points, DimensionalData.Dimensions.Lookups.NoMetadata}}}, Tuple{}, Matrix{Float64}, DimensionalData.NoName, DimensionalData.Dimensions.Lookups.NoMetadata}}, Symbol, Dict{Symbol, Any}}:
ERROR: UndefRefError: access to undefined reference
Stacktrace:
  [1] getindex
    @ ./essentials.jl:13 [inlined]
  [2] getindex
    @ ~/.julia/packages/DimensionalData/BZbYQ/src/array/indexing.jl:28 [inlined]
  [3] first
    @ ./abstractarray.jl:452 [inlined]
  [4] show_after(io::IOContext{…}, mime::MIME{…}, A::DimensionalData.DimGroupByArray{…})
    @ DimensionalData ~/.julia/packages/DimensionalData/BZbYQ/src/groupby.jl:57
  [5] show(io::IOContext{…}, mime::MIME{…}, A::DimensionalData.DimGroupByArray{…})
    @ DimensionalData ~/.julia/packages/DimensionalData/BZbYQ/src/array/show.jl:17
  [6] (::REPL.var"#55#56"{REPL.REPLDisplay{REPL.LineEditREPL}, MIME{Symbol("text/plain")}, Base.RefValue{Any}})(io::Any)
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:273
  [7] with_repl_linfo(f::Any, repl::REPL.LineEditREPL)
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:569
  [8] display(d::REPL.REPLDisplay, mime::MIME{Symbol("text/plain")}, x::Any)
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:259
  [9] display
    @ ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:278 [inlined]
 [10] display(x::Any)
    @ Base.Multimedia ./multimedia.jl:340
 [11] #invokelatest#2
    @ ./essentials.jl:892 [inlined]
 [12] invokelatest
    @ ./essentials.jl:889 [inlined]
 [13] print_response(errio::IO, response::Any, show_value::Bool, have_color::Bool, specialdisplay::Union{…})
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:315
 [14] (::REPL.var"#57#58"{REPL.LineEditREPL, Pair{Any, Bool}, Bool, Bool})(io::Any)
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:284
 [15] with_repl_linfo(f::Any, repl::REPL.LineEditREPL)
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:569
 [16] print_response(repl::REPL.AbstractREPL, response::Any, show_value::Bool, have_color::Bool)
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:282
 [17] (::REPL.var"#do_respond#80"{…})(s::REPL.LineEdit.MIState, buf::Any, ok::Bool)
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:911
 [18] #invokelatest#2
    @ ./essentials.jl:892 [inlined]
 [19] invokelatest
    @ ./essentials.jl:889 [inlined]
 [20] run_interface(terminal::REPL.Terminals.TextTerminal, m::REPL.LineEdit.ModalInterface, s::REPL.LineEdit.MIState)
    @ REPL.LineEdit ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/LineEdit.jl:2656
 [21] run_frontend(repl::REPL.LineEditREPL, backend::REPL.REPLBackendRef)
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:1312
 [22] (::REPL.var"#62#68"{REPL.LineEditREPL, REPL.REPLBackendRef})()
    @ REPL ~/.julia/juliaup/julia-1.10.4+0.x64.linux.gnu/share/julia/stdlib/v1.10/REPL/src/REPL.jl:386
Some type information was truncated. Use `show(err)` to see complete types.

```

---

<div class="post-metadata">

**Author:** ![Raf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raf/32/3383_2.png) [@Raf](https://discourse.julialang.org/u/Raf)\
**Post date:** [August 6, 2024, 3:27pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/32 "2024-08-06T15:27:56Z")

</div>

Yeah that’s just `show` failing with the nested arrays. You can put a semicolon after the line to avoid it.

Maybe make a GitHub issue with your mwe and we can fix that too.

---

<div class="post-metadata">

**Author:** ![aasdelat](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/aasdelat/32/45268_2.png) [@aasdelat](https://discourse.julialang.org/u/aasdelat)\
**Post date:** [August 6, 2024, 3:30pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/33 "2024-08-06T15:30:40Z")

</div>

Ok, I will open the issue.

Additionally, I have observed that, for the example given above, the non threaded solution is faster than the threaded one, but for my data, the threaded solution is faster.

---

<div class="post-metadata">

**Author:** ![Raf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raf/32/3383_2.png) [@Raf](https://discourse.julialang.org/u/Raf)\
**Post date:** [August 6, 2024, 4:04pm UTC](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607/34 "2024-08-06T16:04:30Z")

</div>

Yeah, that’s just threading generally.

I also tested on a much larger array to check it was worthwhile.

[Previous page](https://discourse.julialang.org/t/multithreaded-broadcast-on-a-grouped-dimensionaldata-array/117607.md?page=1)