[Feature Request]: Box plots

[rickpeltier]

Wonder how difficult it would be implement boxplots on maps, similar to what you already produce, but limited to polarPlots. For example, one might have a network of monitors across a region, and may wish to plot annual mean concentrations by year, and/or perhaps an all-year site mean concentration; and do this at several sites across a city. For that matter, I suppose one could call any of the openair or ggplot functions like timeVariation, timePlot, etc.

[jack-davison]

Hi Rick,

You're right - while we currently only support the seven openair directional analysis plots (polar plots, polar annuli, wind roses, etc.) we could in principle use any plot type. All openairmaps does "under the hood" is save the plot images to a temporary directory and then use them as marker images.

In fact, I have a blog post outlining how to do this exact thing:

https://jack-davison.github.io/posts/2023-09-26-ILL-PlotMarkers/

The current suite of plots lend themselves to being used as map markers because they are directional in nature - so allow for things like triangulating sources. But it could be on the cards to have other plot types like boxplots to facilitate easy comparisons between measured concentrations at different sensors.

I think my preference would be for relatively simple graphics (like those of https://github.com/rte-antares-rpackage/leaflet.minicharts) - using a whole timeVariation plot as a marker may end up looking too busy.

Will give this some thought and happy to hear any ideas!

[jack-davison]

Example mock-up below. It's very possible. Some notes:

• I'd be weary about going much more complex than a boxplot, some ideas could be: boxplots, bar charts (e.g., showing mean/max/etc.), simple trend lines (e.g., an annual trend).
• Approaching from an angle of "which other {openair} functions could be used here?" - can perhaps do a simplified TheilSen(), calendarPlot() (switch between months with layer control?), timePlot() and perhaps a single panel of timeVariation() (switch between panels with layer control?).
• The boxplots in the examples below have outliers removed - that's the sort of thing that could be defined by the user in the function args.
• The boxplots below have different pollutants on the y-axis - this could be different (e.g., a single pollutant over years/months/any {openair} "type" arg).
• The boxplots below share a y-axis - could all have distinct y axes? Option for users? Could normalise in some way for multiple pollutants (e.g., so all PM10 concs are cast to [0,1])?
• Basemaps are going to be limited - Carto maps w/out labels seem to be the most effective.

library(openairmaps)
library(openair)
library(leaflet)
library(tidyverse)

# arbitrary postcode
latlng <- convertPostcode("SW1A 1AA")

# get codes near postcode
codes <-
  searchNetwork(
    lat = latlng$lat,
    lng = latlng$lng,
    n = 10,
    map = FALSE,
    year = 2018
  )

# import data
aq <- importUKAQ(codes$code, year = 2018, meta = TRUE)

# pivot longer
aq_long <-
  aq |>
  select(site, code, date, no2, pm10, pm2.5, o3, latitude, longitude) |>
  pivot_longer(no2:o3)

# nest
aq_nested <-
  nest_by(aq_long, site, code, latitude, longitude)

# get a colour palette
vars <- unique(aq_long$name)
cols <- setNames(openair::openColours(scheme = "Spectral", n = length(vars)), vars)

# work out the ymax of boxplots w/out outliers
gg <-
  aq_long |>
  ggplot(aes(x = name, y = value, color = site)) +
  geom_boxplot(outliers = FALSE)

ymax <- max(pretty(ggplot_build(gg)$data[[1]]$ymax))

# function to construct boxplots
make_boxplot <- function(data, ymax) {
  data |>
    ggplot(aes(x = name, y = value)) +
    geom_crossbar(
      data = distinct(data, name),
      aes(
        xmin = name,
        xmax = name,
        ymin = 0,
        ymax = 0,
        y = 0
      ),
      linewidth = 0.1
    ) +
    geom_boxplot(na.rm = TRUE,
                 aes(fill = name, color = name),
                 alpha = 0.5, 
                 outliers = FALSE,
                 linewidth = 0.5) +
    theme_void() +
    theme(legend.position = "none", aspect.ratio = 1) +
    scale_color_manual(values = cols, aesthetics = c("fill", "color")) +
    scale_y_continuous(limits = c(0, ymax))
}

# temporary directory to save plots
t <- tempdir()

# create plots & urls
plots <-
  aq_nested |>
  mutate(plot = list(make_boxplot(data, ymax)),
         url = paste0(t, "/", latitude, "_", longitude, ".png"))

# save plots
purrr::walk2(.x = plots$plot,
             .y = plots$url,
             .progress = TRUE,
             .f = ~ {
               ggsave(
                 filename = .y,
                 plot = .x,
                 width = 1,
                 height = 1,
                 dpi = 300,
                 bg = "transparent"
               )
             })

# make leaflet map
leaflet(plots) |>
  addProviderTiles(providers$CartoDB.PositronNoLabels) |>
  addMarkers(icon = ~ makeIcon(url, iconWidth = 75, iconHeight = 75),
             label = ~ site) |>
  addLegend(colors = cols,
            labels = quickTextHTML(names(cols)),
            title = "Pollutant")

[rickpeltier]

yes, this is just what I was asking about! And I fully recognize that this could become really complicated really quickly if you used complex data viz outputs as your markers (like the full timeVariation() output - which would be a true chartjunk). But we often use these outputs as infographics to disseminate data to a community, so being able to do so natively would be really powerful. Thanks!

[jack-davison]

I think we're on the same page, then!

I've just released a new version of {openairmaps} with some updates I needed to do for static map generation, so it's timely to think about "what's next?" - I'll have a good think about how to implement this because I think it will be useful - especially for things like sensor networks where there are lots of measurement locations.