Removing white border plotting maps with geom_sf

I am attempting to remove the white borders around a facet plot that I created with geom_sf(). Previously I was using the now deprecated rgdal package and coord_cartesian(). I am trying to get a similar graph plotted using geom_sf() and coord_sf(), but can't seem to get it to work.

My code:

land <- sf::read_sf("./shapefiles/Grey_WA_OR.shp")
land <- sf::st_transform(land, crs = "+proj=longlat +datum=WGS84")

water <- sf::read_sf("./shapefiles/NW_Marine.shp")
water <- sf::st_transform(water, crs = "+proj=longlat +datum=WGS84")

heatmap <- ggplot()+ theme_bw()+
  geom_sf(data = water, fill = "paleturquoise3", color = "turquoise4") +
  geom_sf(data = land, fill='grey96') +
  coord_sf(xlim = c(-126.3, -122.3), ylim = c(48.3, 50.39)) +
  geom_point(data= CPUE.bymonth.mean.nonzero, mapping = aes(x=Lon, y=Lat, 
            size=MeanCPUE), alpha = 0.5, colour= "coral", stroke = 0.7)+  
  geom_point(data = CPUE.bymonth.mean.zeros, mapping = aes(x = Lon, y = Lat, 
                          size = MeanCPUE), alpha = 0.7, color = "gray25") +
  scale_size_continuous(breaks = c(0, 500,1000, 3000), name = "Mean CPUE", 
                        range = c(2.5,15)) +
  theme(axis.text.x = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks = element_blank(),
        axis.title.x = element_blank(), 
        axis.title.y = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank())+
  guides(size = guide_legend(override.aes = list(color = 
        c("gray25","coral","coral", "coral"), alpha = 0.6), size = 12))+
  facet_grid(Year ~ Month)+
  theme(
    strip.text.x = element_text(size = 12, color = "black", face = "bold"),
    strip.text.y = element_text(size = 12, color = "black", face = "bold"))

ggsave(filename = "figures/CPUE_map.png", plot = heatmap, 
       width = 12, height = 6, dpi = 600 )

What I get:

What I hope to get:

I've looked at the cowplot package but haven't been able to get that to work unfortunately, not sure if that package holds the solution? I've also tried removing the margins with theme(plot.margin = unit(c(0,0,0,0), "cm")), but no luck.

The issue lies with you maintaining the same width and height of your first ggsave() output even though you have projected your data to WGS84/EPSG:4326. Note how Vancouver Island is much more elongated in the first image relative to the second image. That's because coord_sf() automatically plots data to the relative coordinates of an sf object. With coord_cartesian(), the x and y values are treated as equal. If you need to fill the plot area, you have three options:

1. Maintain the same xlim/ylim values in coord_sf() and simply adjust the width and height values in ggave()
2. Edit the xlim/ylim values in coord_sf() until the output 'fits' to the existing width and height values
3. Convert your sf objects to dataframes and plot them using geom_polygon() in combination with coord_fixed()

This simplified reprex demonstrates option 1 and 3, the latter I believe is the one you want. As you did not provide any sample data, I have used the ne_countries() dataset from the rnaturalearth package in place of "land". To keep it simple, I only created sample points to represent "CPUE.bymonth.mean.nonzero" and used panel.background = in place of "water".

Option 1

library(rnaturalearth)
library(sf)
library(dplyr)
library(ggplot2)

# Load world map data
land <- ne_countries(scale = "medium")

# Create sample point data
set.seed(1)
CPUE.bymonth.mean.nonzero <- data.frame(Lon = runif(10, -126, -124),
                                        Lat = runif(10, 48.5, 49.5),
                                        Year = rep(c(2022, 2023), each = 5),
                                        Month = rep(month.abb[4:8], 2),
                                        MeanCPUE = sample(100:4000, 10))

# Plot
heatmap <- ggplot() +
  geom_sf(data = land) +
  geom_point(data= CPUE.bymonth.mean.nonzero,
             aes(x = Lon, y = Lat, size = MeanCPUE),
             alpha = 0.5, 
             colour = "coral",
             stroke = 0.7) +
  coord_sf(xlim = c(-126.2986, -122.2988), 
           ylim = c(48.30015, 50.39008)) +
  scale_size_continuous(breaks = c(0, 500,1000, 3000), 
                        name = "Mean CPUE", 
                        range = c(2.5, 15)) +
  theme_bw() +
  theme(axis.text.x = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks = element_blank(),
        axis.title.x = element_blank(), 
        axis.title.y = element_blank(),
        panel.background = element_rect(fill = "paleturquoise3", colour = NA),
        legend.key = element_rect(fill = "white", colour = NA),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank()) +
  guides(size = guide_legend(override.aes = list(color = c("coral", "coral", "coral"), 
                                                 alpha = 0.6), size = 12)) +
  facet_grid(Year ~ Month) +
  theme(strip.text.x = element_text(size = 12, color = "black", face = "bold"),
        strip.text.y = element_text(size = 12, color = "black", face = "bold"))

ggsave("figures/CPUE_map.png",
       heatmap, 
       width = 12,
       height = 3.8,
       dpi = 600)

Here's the result, note that this image is only 150dpi:

Option 3

In this example, I have used coord_fixed(ratio = 2.5), which seems to be a close approximation of your original plot. Adjust if necessary.

# Convert land sf to df
land_df <- land |>
  select(geometry) |>
  st_cast("POLYGON") |>
  st_coordinates() |>
  as.data.frame() |>
  rename(lon = X, lat = Y)

# Plot
heatmap <- ggplot() +
  geom_polygon(data = land_df, 
               aes(x = lon, y = lat, group = L2),
               fill = "grey96",
               colour = "turquoise4") +
  geom_point(data= CPUE.bymonth.mean.nonzero,
             aes(x = Lon, y = Lat, size = MeanCPUE),
             alpha = 0.5, 
             colour = "coral",
             stroke = 0.7) +
  coord_fixed(xlim = c(-126.2986, -122.2988), 
              ylim = c(48.30015, 50.39008),
              ratio = 2.5) +
  scale_size_continuous(breaks = c(0, 500,1000, 3000), 
                        name = "Mean CPUE", 
                        range = c(2.5, 15)) +
  theme_bw() +
  theme(axis.text.x = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks = element_blank(),
        axis.title.x = element_blank(), 
        axis.title.y = element_blank(),
        panel.background = element_rect(fill = "paleturquoise3", colour = NA),
        legend.key = element_rect(fill = "white", colour = NA),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank()) +
  guides(size = guide_legend(override.aes = list(color = c("coral", "coral", "coral"), 
                                                 alpha = 0.6), size = 12)) +
  facet_grid(Year ~ Month) +
  theme(strip.text.x = element_text(size = 12, color = "black", face = "bold"),
        strip.text.y = element_text(size = 12, color = "black", face = "bold"))

ggsave("figures/CPUE_map.png",
       heatmap, 
       width = 12,
       height = 6,
       dpi = 600)

Again, this example plotted at 150dpi, but the width and height values in ggsave() match your original plot.