Utilizing https://github.com/geoarrow/deck.gl-layers for mapgl

[tim-salabim]

Hi @walkerke and everyone who may be interested (and tagged :-) ).

As alluded to in #57 (comment) I want to implement something that utilizes the blazing fast data transfer from R to the browser using apache arrow. With the awesome work of @paleolimbot (nano-/geoarrow) & @kylebarron (geoarrow/deck.gl-layers) we can leverage this nicely now from R.

My earlier efforts trying to add this on top of leaflet are not too encouraging, so I have been playing around with the maplibre part of mapgl and things look quite promising. Below is a pseudo-reproducible example (reading the data directly from the github source takes awefully long, hence the GDAL download part to get the data) that should highlight what this is capable of. The roughly 1M building footprints of Idaho, US render nicely and map creation is basically instant!

A few things I noticed:

• map handling is not quite as pleasant as when using geojson layers (as they seem to be tiled on-the-fly by maplibre), but that is probably something that can be implemented in deckgl in the future (@felixpalmer) - at the cost of speed at which the map is ready?
• globe view rendering in deckgl is not quite stable yet (try zooming in and out when switched to globe view using the button in the top-right corner)
• I would've set the zoom and bounds appropriately (to Idaho in the example below) using set_view but it seems it only respects the longitude part of the center and lat is always 0
• the example below includes a tiny script to add a globeControl button, which should probably be in mapgl itself (I can prepare a PR if desired).

There are probably a lot of things left unanswered here, but I wanted to get the ball rolling and start discussion on how to proceed from here. One obvious question is whether this should eventually make its way into mapgl's code base or if its better hosted in a separate package. I clearly prefer the latter!

All in all, I am quite happy with how this is working. Feedback and suggestions more than welcome!

library(mapgl)
library(geoarrow)
library(nanoarrow)
library(sf)

source("https://gist.githubusercontent.com/tim-salabim/20f330a442f4861a1620711bc04bf41d/raw/6f21f0e29be14b1a6ca439ca08b42ecf1f6f54ee/add_deckgl_polygons.R")

## Download the building footprint data using ogr2ogr in your terminal:
## see overview of footprints here:
## https://github.com/microsoft/USBuildingFootprints

# ogr2ogr \
# ~/Downloads/data.gpkg \
# /vsizip//vsicurl/https://minedbuildings.z5.web.core.windows.net/legacy/usbuildings-v2/Idaho.geojson.zip

data = st_read("~/Downloads/data.gpkg")

## generate map 
m = maplibre(style = 'https://basemaps.cartocdn.com/gl/positron-gl-style/style.json') |>
  add_navigation_control(visualize_pitch = TRUE) 

## now add the dependencies - couldn't find a way to do this directly loading from gists
## you probably need to adjust the file.path if not on Linux
jsDir = file.path("~/Downloads", "jsDeps")
dir.create(jsDir)
globeJsFile = file.path(jsDir, "addGlobeControl.js")
deckglPolygonsFile = file.path(jsDir, "addDeckglPolygons.js")

download.file(
  "https://gist.github.com/tim-salabim/2db92c7525b7623222aadf20bd7f4dcf/raw"
  , destfile = globeJsFile
)

download.file(
  "https://gist.github.com/tim-salabim/8371fa6523c689bbd3c597a1f6334dfe/raw"
  , destfile = deckglPolygonsFile
)

m$dependencies = c(
  m$dependencies
  , list(
    htmltools::htmlDependency(
      name = "globeControl"
      , version = "0.0.1"
      , src = jsDir
      , script = "addGlobeControl.js"
    )
  )
  , list(
    htmltools::htmlDependency(
      name = "deckglPolygons"
      , version = "0.0.1"
      , src = jsDir
      , script = "addDeckglPolygons.js"
    )
  )
)
  
## my RStudio viewer has had issues for a while being terribly slow, so view in browser
options(viewer = NULL)

m |> 
  add_deckgl_polygons(
    data = data 
    , layerId = "test"
  ) |>
  add_layers_control(collapsible = TRUE, layers = c("test"))

[cboettig]

deck.gl support would be cool. A side comment, but PMTiles works rather nicely for plotting millions of polygons, like this building dataset. This should be a reproducible example that runs quickly:

library(mapgl)

pmtiles <-"https://data.source.coop/cholmes/overture/overture-buildings.pmtiles"
m <-
  maplibre(center=c(-74.0, 40.7), zoom=14, pitch=60, bearing=-17) |>
  add_vector_source("source",
                    url = paste0("pmtiles://", pmtiles)
  ) |>
  add_fill_extrusion_layer(
    id = "layer",
    source = "source",
    source_layer = "buildings",
    fill_extrusion_color = interpolate(
      column = "height",
      values = c(0, 200, 400),
      stops = c("lightgray",  "royalblue", "lightblue")
    ),
    fill_extrusion_height = list("get", "height"),
    tooltip = "height"
  )

m

should have buildings of the world.

[tim-salabim]

Yes! That is awesome.
It's a different use case though IMO. Vector tiles are really the only option for larger than memory data, but take time to prepare... The deckgl approach is for, let's say, medium sized data that still fits in memory, but is too large for reasonably being handled by a geojson approach.

[paleolimbot]

This is the dream! Let me know if anything comes up from the geo/nanoarrow end of things that would be helpful for me to take a look at!

[paleolimbot]

Ooof, that's a big one (but you're right that it's super useful!). I'm shooting for a nanoarrow release in Feb/March and will try to get there!

[tim-salabim]

Well, as I said, it's not a biggie for my side as converting factor to character is easy on the R side. So no worries if it's not a priority

[tim-salabim]

@paleolimbot also a solution to this would be great, as we don't need the overhead of coordinates that are too precise for the webmap.

[paleolimbot]

It's definitely a bug and should be fixed, although I am not sure it will help reduce the size of the IPC you need to send. Rounding coordinates does tend to improve the compression ratio you're able to achieve (but I am not sure that will improve speed...in my brief experience the bandwidth to localhost is fast enough that compression doesn't help). If we could use floats that might help (but if we added that to geoarrow we'd have an absolute gongshow of memory layouts to deal with).

[tim-salabim]

Thanks for the quick reply. Even if it doesn't help data transfer performance, my hope is that it will help rendering performance of the deckgl layer on the browser side. In any case, the data I am playing with ends up with 15 digit geographic coordinates (cadastral data) and we all know that this is a bit over the top

https://xkcd.com/2170/

[cboettig]

💯 agree with @tim-salabim -- vector tiles are ideal for this very big data we can take the time to pre-render and host, but not a viable option for data assets that might be generated on the fly, e.g. in a shiny app. (aside but I also ❤️ that toggle to globe projection you have there!)

I would definitely love to see deck-gl layers! And being able to leverage arrow would be amazing.

Meanwhile, there are some decent if largely intuitive suggestions in the maplibre docs, as well as extensions that support other formats than geojson.

FWIW, I've been spitballing how to do quick render and aggregates for very large data using @walkerke 's amazing maplibre bindings here and a few extensions. Here's a rough proof-of-concept that aggregates species occurrences to arbitrary zoom levels (h3 hex sizes 0-11) over arbitrary selected spatial areas from a global database of 3 billion rows (347 GB) using duckdb + h3 extension + mapgl https://shiny.carlboettiger.info/

[tim-salabim]

This is cool! My dream is to have something like this automatically in the browser (with some pre-defined arguments defined in R probably). I think this is all doable now, given that we can send all these data to the browser pretty much instantly. Having a switch between an aggregated layer like H3 and the original points layer should be doable without too much hassle, I guess. A good and concise API may be a challenge though. The harder issue are line and polygon layers, as we need some clear heuristics on how to aggregate (more likely thin out) the data depending on zoom.
Yet, we could start by trying something like "only show the largest x polygons within view bbox" or something similar. I think there are solutions to do calculations like this in the browser (duckdb-wasm for sure, but there are other libs as well IIRC). Just thinking out loud here...

[walkerke]

This is amazing stuff - thanks for the discussion @tim-salabim @cboettig @paleolimbot!

Leveraging arrow + deck.gl in mapgl would be fantastic, especially for these big data visualization use-cases. Mapbox / MapLibre tend to choke once you hit a million or so GeoJSON features.

With respect to whether this merits a separate package, or should be included in mapgl's codebase: I'll give this some thought, I'm open to either. If it were to be included, we'd probably want some degree of consistency with respect to how people interact with and style regular Mapbox / MapLibre layers vs. deck.gl layers, which would take some work. That said, mapgl has a fairly well-developed Shiny interface that I could see people wanting to use with the deck.gl layers.

[tim-salabim]

If it were to be included, we'd probably want some degree of consistency with respect to how people interact with and style regular Mapbox / MapLibre layers vs. deck.gl layers, which would take some work.

I would suggest to stay as close to deckgl specifications/lingo as possible. It just makes more sense to me to mirror all arguments etc from the original implementation, so that it is easy to look up/link to respective online documentation (which is fairly good for deckgl).

That said, mapgl has a fairly well-developed Shiny interface that I could see people wanting to use with the deck.gl layers.

Yeah, shiny is always a layer on top of the original, so there's room for diverting from the original definitions. But for adding a deckgl suite of add* functions, I guess it is mostly internal stuff that needs to be added (I could be wrong here)?

My suggestion is to start with a package, and we can take it from there. How does that sound?

[tim-salabim]

Another point to consider is that utilising geoarrow is dependency heavy, both on the R ans JS side

[tim-salabim]

Package now live at https://github.com/tim-salabim/geoarrowDeckgl
So far, we have Scatterplot and Polyfgon layers (first drafts). Others to follow as soon as those are fully functional.

[walkerke]

Another point to consider is that utilising geoarrow is dependency heavy, both on the R ans JS side

Makes sense to me!

[paleolimbot]

Perhaps worth considering rewriting geoarrow-R in C to keep the dependency cost down! (/sarcasm...is it actually an issue in the presence of something like sf?)

[tim-salabim]

is it actually an issue in the presence of something like sf?

I don't think so. It's benefits are well worth the dependency heaviness. At least IMHO

[cboettig]

Love the vision @tim-salabim outlines here! Re:

The harder issue are line and polygon layers, as we need some clear heuristics on how to aggregate (more likely thin out) the data depending on zoom.

Just a note on this, some approaches to these heuristics are already outlined in the tippecanoe library (C++), which builds vector tiles like PMTiles. For instance, line-type data can be filtered by weight at different zooms (i.e. low zooms show highways, only high zooms include tiny lanes and streets). With polygons, the strategies involve simplifying polygons at low zooms, and dropping very small polygons. The right choice depends on the data specifics though -- e.g. dropping polygons is good for the building footprint example above but looks weird if the polygons should tessellate the area (like census tracts, say).

(total tangent but R bindings to tippecannoe would be awesome).

Cool to see geoarrowDeckgl!

[tim-salabim]

(total tangent but R bindings to tippecannoe would be awesome).

I agree that having R bindings to tippecanoe would be awesome, but just to be clear, I was thinking about doing this on-the-fly in the browser as the user zooms and pans around...

[kylebarron]

I've thought about on-the-fly server-generated (i.e. from the Python process) tiles using https://github.com/mapbox/supercluster (i.e. https://github.com/developmentseed/supercluster-rs) or https://github.com/mapbox/geojson-vt (to be ported to something based on Arrow). But haven't ever truly tested it.

[tim-salabim]

I think geoarrow-vt would be super awesome!

[tim-salabim]

And while I think marker/point clustering like in supercluster is useful, a more general approach IMO is to cluster onto a fixed grid/hexlayer with a user-definable aggregation function (e.g. count, max, min, mean or whatever is best suited for the data and analysis at hand).

[felixpalmer]

@tim-salabim take a look at https://loaders.gl/examples/tiles/table-tiler I think it does exactly what you are after

[tim-salabim]

Thanks Felix! Yes. that looks like what we want. How would we use this with geoarrow data instead of geojson, though?

[kylebarron]

I think it's just a wrapper around geojson-vt. We'd need to rewrite/update it to use Arrow under the hood.

[tim-salabim]

That's what I thought, thanks

[felixpalmer]

Yes, you would need to update it to support Arrow as @kylebarron said