Search code examples
rggplot2dplyrmappinglatitude-longitude

How to plot a map based on actual latitude and longitude?


I have a region divided by longitude and latitude. I know that the distance of 1 degree of latitude is 111,320m, and the distance of 1 degree of longitude varies with the latitude. The size of my grid is 25m×25m, and the grid shown in the figure (a total of 600 grids) is too large (see the figure below), which is not realistic.

enter image description here

My code is as below:

library(ggplot2)
library(dplyr)
library(readr)
color_breaks <- c(0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0)
colors <- c("#E31A1C", "#FF7F00", "#FDBF6F", "#E9E4A6", "#A4D4A9", "#B2DF8A", "#33A02C", "#1F78B4")

p1 <- ggplot(merged_data, aes(x = longitude, y = latitude, z = median_mafruit)) +
  labs(
    title = "Data1",
    subtitle = "10 years",
    x = "Longitude (°)",
    y = "Latitude (°)",
    fill = expression("Amount")
  ) +
  theme_minimal()

p1<-p1 +
  stat_summary_2d(
    aes(fill = after_stat(
      cut(
        value,
        breaks = c(0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0),
        include.lowest = TRUE,
        right = FALSE
      )
    )),
    bins = 30,
    show.legend = TRUE
  ) +
  scale_fill_manual(
    values = colors,
    drop = FALSE,
    guide = guide_legend(reverse = TRUE)
  ) +
  theme(
    plot.title = element_text(size = 14, face = "bold"),
    plot.subtitle = element_text(size = 12, face = "italic"),
    legend.title = element_text(size = 10),
    legend.text = element_text(size = 8),
    legend.key.size = unit(0.5, "cm"),
    legend.spacing.x = unit(0.5, "cm"),
    legend.spacing.y = unit(0.5, "cm")
  )
p1

Part of my dataset is shown below:

structure(list(Order = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 
13, 14, 15, 16, 17, 18, 19, 20), latitude = c(43.2143, 43.3697, 
43.3909, 43.3926, 43.3961, 43.3978, 43.066, 43.2215, 43.368, 
43.435, 43.4434, 43.4623, 43.4895, 43.4856, 43.3738, 43.4761, 
43.5102, 43.5118, 43.5062, 43.4933), longitude = c(-4.479, -4.4804, 
-4.4606, -4.4484, -4.4241, -4.412, -4.1046, -4.1049, -4.1031, 
-4.0818, -4.021, -3.9502, -3.8184, -3.7798, -3.7279, -3.7147, 
-3.5971, -3.5849, -3.5585, -3.5177), median_mafruit = c(2.73, 
1.095, 1.115, 2.73, 0.527, 0.527, 0.962, 1.039, 1.039, 2.73, 
2.73, 2.73, 2.73, 2.73, 0.544, 2.73, 2.73, 2.73, 0.478, 2.73)), row.names = c(NA, 
-20L), spec = structure(list(cols = list(Order = structure(list(), class = c("collector_double", 
"collector")), latitude = structure(list(), class = c("collector_double", 
"collector")), longitude = structure(list(), class = c("collector_double", 
"collector")), median_mafruit = structure(list(), class = c("collector_double", 
"collector"))), default = structure(list(), class = c("collector_guess", 
"collector")), delim = ","), class = "col_spec"), class = c("spec_tbl_df", 
"tbl_df", "tbl", "data.frame"))

I want to keep appropriate gaps between grids, and the lengths of the horizontal and vertical coordinates of the image should be more in line with the actual proportion (currently 2 latitudes correspond to 5 longitudes, which is not practical).


Solution

  • You can use coord_sf(), which among other things interprets coordinates relative to a geographic or projected coordinate system. In your example, I have assumed your coordinates are WGS84/EPSG:4326. If they are in a different coordinate system, you can look up the relevant EPSG code online, then replace 4326 with the correct code.

    If you intend to use ggplot2 to plot spatial data in the future, I recommend you explore the sf and terra packages. For sf objects, you can use geom_sf(), which 'intuitively' handles coordinates etc. For SpatRaster objects, you will need the tidyterra package to use geom_spatraster().

    library(ggplot2)
    library(dplyr)
    
    p1 <- ggplot(merged_data, aes(x = longitude, y = latitude, z = median_mafruit)) +
      labs(
        title = "Data1",
        subtitle = "10 years",
        x = "Longitude (°)",
        y = "Latitude (°)",
        fill = expression("Amount")
      ) +
      theme_minimal()
    
    p1 <- p1 +
      stat_summary_2d(
        aes(fill = after_stat(
          cut(
            value,
            breaks = c(0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0),
            include.lowest = TRUE,
            right = FALSE
          )
        )),
        bins = 30,
        show.legend = TRUE
      ) +
      scale_fill_manual(
        values = colors,
        drop = FALSE,
        guide = guide_legend(reverse = TRUE)
      ) +
      coord_sf(datum = 4326
      ) +
      theme(
        plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 12, face = "italic"),
        legend.title = element_text(size = 10),
        legend.text = element_text(size = 8),
        legend.key.size = unit(0.5, "cm"),
        legend.spacing.x = unit(0.5, "cm"),
        legend.spacing.y = unit(0.5, "cm")
      )
    p1
    

    1