I am trying to plot some home ranges, and I'd like to "zoom" to a smaller area. I can plot the range polygon without issue, but the problem occurs when I try to look closer at a specific area.
When I set the X and Y limits to the area of interest (using coord_sf()), the polygon is removed. I think this is because ggplot automatically removes data outside of the plotting area. For regular data, I know coord_cartesian() solves this problem, but that seems not to be an option when using polygons and geom_sf() ("Error: geom_sf() must be used with coord_sf()").
Is there some equivalent for this type of data, or perhaps a simple way to crop the polygon based on X and Y limits?
Example below.
library(sf)
library(ggmap)
library(adehabitatHR)
library(rnaturalearth)
library(lwgeom)
data <- data.frame(x = c(-50.3, -49.9, -50.0, -50.6, -55.3, -55.4, -55.5, -55.3, -54.9, -54.4, -51.5, -51.2, -50.8, -50.3),
y = c(50.3, 48.8, 48.1, 47.4, 48.2, 47.4, 50.1, 48.1, 47.5, 50.7,50.4, 50.7, 50.5, 48.3))
data$id <- as.factor("a")
#create a SpatialPointsDataFrame by defining coordinates
coordinates(data) <- c("x", "y")
# currently there is no CRS for this data. Since it is lat/long, set it as such:
proj4string(data) = CRS("+proj=longlat +datum=WGS84 +no_defs")
# now it needs to be converted to UTM
data <- spTransform(data, CRS("+init=epsg:32621 +proj=utm +zone=21 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 "))
# NOTE: I know proj4strings are no longer current but this is old code that still works for the purpose of this example
# run code to determine minimum convex polygon
mcp <- mcp(data, percent = 95, unin = c("m"),
unout = c("km2"))
# view resulting polygons by plotting them
plot(data, col = as.factor(data@data$id), pch = 16, axes = TRUE)
plot(mcp, col = alpha(1:7, 0.5),add = TRUE)
# now plot polygons on map
# import map and check projection
canada <- ne_countries(country = "canada", scale = "large", returnclass = "sf") # in rnaturalearth package
st_crs(canada)
# convert to matching utm
nl <- st_transform_proj(canada, "+init=epsg:32621 +proj=utm +zone=21 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0")
# now plot larger area
ggplot(data = nl) +
geom_sf(data = nl, fill = "gray80") +
coord_sf(xlim = c(100000, 1200000), ylim = c(5000000, 6100000)) + # map coordinates in utm
labs(x = "Longitude", y = "Latitude", colour = "") + # labels
geom_polygon(data = mcp, fill = "goldenrod", alpha=0.7, aes(x = long, y = lat, group = group))
# zoom to area of interest
ggplot(data = nl) +
geom_sf(data = nl, fill = "gray80") +
coord_sf(xlim = c(600000, 920000), ylim = c(5100000, 5400000)) + # map coordinates in utm
labs(x = "Longitude", y = "Latitude", colour = "") + #labels
geom_polygon(data = mcp, fill = "goldenrod", alpha = 0.7, aes(x = long, y = lat, group = group))
# polygon removed
Best Answer
Remove the alpha argument and you won't see it disappear. I'm not sure why, but I'm not seeing this resolved elsewhere. It also seems to do with what kind of machine you're working on (i.e., it happens on Windows but not Mac OS).