Search code examples
rggplot2spatialr-rastertidyterra

ggplot of raster maps: Horizontal / vertical grid cells distortion


I am plotting multiple projected raster objects (SpatRaster or raster stack; equal area projection; EPSG:3035), all having the same extent. Example map

However, the mapped grid cells are distorted, either horizontally or vertically (i.e. not square), or shown larger than they should be.

# vertically distorted  
ggplot2::ggplot() +
    tidyterra::geom_spatraster(data = terra::rast(ExampleMap[[1]])) +
    paletteer::scale_fill_paletteer_c(na.value = NA, "viridis::plasma")

# horizontally distorted  
ggplot2::ggplot() +
    tidyterra::geom_spatraster(data = terra::rast(ExampleMap[[2]])) +
    paletteer::scale_fill_paletteer_c(na.value = NA, "viridis::plasma")


Here is a full script for plotting one layer using different geoms:

require(raster)
require(stara)
require(ggplot2)

ExampleMap2 <- ExampleMap[[1]] %>% trim()

Using geom_stars

ggplot2::ggplot() +
  stars::geom_stars(data = stars::st_as_stars(ExampleMap2)) + 
  ggplot2::coord_equal() +
  ggplot2::theme_void() +
  ggplot2::scale_fill_continuous(na.value = NA) + 
  guides(fill = guide_legend(NULL)) +
  theme(aspect.ratio = 1) + 
  labs(title = "stars::geom_stars")

Using tidyterra::geom_spatraster

ggplot2::ggplot() +
  tidyterra::geom_spatraster(data = terra::rast(ExampleMap2)) + 
  ggplot2::theme_void() + 
  ggplot2::scale_fill_continuous(na.value = NA) + 
  guides(fill = guide_legend(NULL)) +
  theme(aspect.ratio = 1) + 
  labs(title = "tidyterra::geom_spatraster")

Using ggplot2::geom_raster

DT <- as.data.frame(ExampleMap2, xy = TRUE)

ggplot2::ggplot() + 
  ggplot2::geom_raster(aes(x = x, y = y, fill = as.numeric(as.character(ExampleMap1))), data = DT) + 
  ggplot2::theme_void() +
  ggplot2::scale_fill_continuous(na.value = NA) + 
  guides(fill = guide_legend(NULL)) + 
  coord_equal() +
  theme(aspect.ratio = 1) + 
  labs(title = "ggplot2::geom_raster")

Using ggplot2::geom_tile [this seems to work as expected]

ggplot2::ggplot() + 
  ggplot2::geom_tile(aes(x = x, y = y, fill = as.numeric(as.character(ExampleMap1))), data = DT) + 
  ggplot2::theme_void() +
  ggplot2::scale_fill_continuous(na.value = NA) + 
  guides(fill = guide_legend(NULL)) + 
  coord_equal() + 
  labs(title = "ggplot2::geom_tile")

What is the reason for this distortion and how to ensure that this does not happen using e.g. tidyterra::geom_spatraster or ggplot2::geom_raster?

Thanks


Solution

  • I modified my original answer, noted now under Old Answer

    It seems to be related with how ggplot handles NA on geom_raster() when there are few values left (see https://github.com/tidyverse/ggplot2/issues/5523 and https://github.com/tidyverse/ggplot2/issues/3025. Instead use na.value = "transparent"

    # Get the data
    
    tmp <- "ExampleMap.zip"
    download.file(
      "https://github.com/dieghernan/tidyterra/files/13322139/ExampleMap.zip",
      tmp
      
    )
    dataset <- unzip(tmp, list = TRUE)[1]
    unzip(tmp)
    
    
    load("ExampleMap.RData")
    
    library(terra)
    #> terra 1.7.55
    library(raster) # For the reprex
    #> Loading required package: sp
    library(tidyterra)
    #> 
    #> Attaching package: 'tidyterra'
    #> The following object is masked from 'package:raster':
    #> 
    #>     select
    #> The following object is masked from 'package:stats':
    #> 
    #>     filter
    library(ggplot2)
    library(paletteer)
    
    ExampleMap2 <- rast(ExampleMap[[2]])
    
    # Delete some extra cells for clarity
    smp <- ExampleMap2 %>%
      slice(380000:400000, .keep_extent = FALSE) %>%
      trim()
    
    
    # Use transparent instead of NA
    
    ggplot() +
      geom_spatraster(data = smp) +
      paletteer::scale_fill_paletteer_c(na.value = "transparent", "viridis::plasma")
    

    Created on 2023-11-15 with reprex v2.0.2

    Old answer

    I added several tests to tidyterra and didn't find an issue with this. I think here the problem is not that tidyterra is distorting the cells, but how the ggplot2 palette is filling the plot, and related with the value of na.value of the fill scale.

    Maybe this can be an issue on ggplot2?

    See how the ggplot2 is affected by this parameter, and a workaround using the geom_tile() approach:

    
    
    # Get the data
    
    tmp <- tempfile(fileext = ".zip")
    download.file(
      "https://github.com/dieghernan/tidyterra/files/13322139/ExampleMap.zip",
      tmp
    )
    dataset <- unzip(tmp, list = TRUE)[1]
    unzip(tmp, exdir = tempdir())
    
    library(raster)
    
    ## Loading required package: sp
    
    load(file.path(tempdir(), dataset))
    
    library(terra)
    
    ## terra 1.7.55
    
    library(tidyterra)
    library(ggplot2)
    library(paletteer)
    
    ExampleMap2 <- rast(ExampleMap[[2]])
    
    # Delete some extra cells for clarity
    as_tibble(ExampleMap2, cells = TRUE, na.rm = TRUE)
    
    ## # A tibble: 7 × 2
    ##     cell ExampleMap2
    ##    <int>       <dbl>
    ## 1 281039           1
    ## 2 343040           1
    ## 3 389362           1
    ## 4 391042           1
    ## 5 392722           1
    ## 6 392725           1
    ## 7 435454           1
    
    smp <- ExampleMap2 %>%
      slice(380000:400000, .keep_extent = FALSE) %>%
      trim()
    
    # With terra
    plot(smp)
    

    enter image description here

    # With tidyterra, no palettes
    
    p <- ggplot() +
      geom_spatraster(data = smp)
    
    # No distortion
    p
    

    enter image description here

    # When I add the palette
    p +
      paletteer::scale_fill_paletteer_c(na.value = NA, "viridis::plasma")
    

    enter image description here

    The problem is the na.value param

    p +
      paletteer::scale_fill_paletteer_c(na.value = "grey80", "viridis::plasma")
    

    enter image description here

    # Full script
    
    ExampleMap1 <- rast(ExampleMap[[1]]) %>% trim()
    p2 <- ggplot() +
      geom_spatraster(data = ExampleMap1, maxcell = Inf)
    
    p2 +
      paletteer::scale_fill_paletteer_c(na.value = NA, "viridis::plasma")
    

    enter image description here

    # But
    p2 +
      paletteer::scale_fill_paletteer_c(na.value = "grey99", "viridis::plasma")
    

    enter image description here

    Workaround

    ExampleMap1 <- rast(ExampleMap[[1]]) %>% trim()
    
    # Make use of the fortify property to pass to geom_tile
    ggplot(ExampleMap1, aes(x,y)) +
        geom_tile(aes(fill=ExampleMap1)) +
    # Palette
        paletteer::scale_fill_paletteer_c(na.value = NA, "viridis::plasma") +
    # Add spatial features to the plot
        coord_sf(crs=terra::crs(ExampleMap1)) +
        labs(x="", y = "")
    
    

    enter image description here