Medición del porcentaje de solapamiento entre los rangos de origen y otro archivo shape en R

Question

Dispongo de datos de seguimiento de varios animales y he creado estimaciones de densidad de núcleos de sus áreas de distribución. Quiero averiguar la proporción de cada área de distribución que se encuentra dentro de zonas protegidas. Utilizo el método adehabitatHR para calcular los rangos domésticos y el wdpar para descargar las zonas protegidas de la Base de Datos Mundial sobre Zonas Protegidas. Pero no estoy seguro de cómo vincular los dos para calcular la proporción de solapamiento. He aquí un ejemplo del código que tengo hasta ahora:

library(amt)
library(wdpar)
library(rgdal)
library(adehabitatHR)

data(deer)
mini_deer <- deer[1:100, ]
mini_deer$id <- ifelse(mini_deer$x_ < 4313000, "A", "B")

mini_deer <-
  mk_track(
    mini_deer,
    .x = x_,
    .y = y_,
    .t = t_,
    id = id,
    crs = CRS(
      "+init=epsg:3035 +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"
    )
)

# Keep only the id column and the x and y coordinates
mydata.sp <- mini_deer[, c("id", "x_", "y_")] 
coordinates(mydata.sp) <- c("x_", "y_")

# kde
kernel.ref <- kernelUD(mydata.sp, h = "href")  # href = the reference bandwidth
mydata.kernel.poly <- getverticeshr(kernel.ref, percent = 95, unout = "m2") 
print(mydata.kernel.poly)  # returns the area of each polygon

# get the protected area data
ger_raw_data <- wdpa_fetch("Germany", wait = TRUE)
pa_data <- st_transform(ger_raw_data, 3035)

# scene missing
# proportion of overlap for A and B

Answer 1

En este caso, parece que no hay solapamiento entre sus hábitats (en rojo) y las zonas protegidas (en naranja). Este mapa se ha elaborado en QGIS:

pero en general una vez que tienes ambos conjuntos de datos como sf objetos de clase y como polígonos (hay algunos puntos en los datos del área protegida) entonces se hace st_intersection del hábitat con la zona protegida para obtener la zona de solapamiento. Para sus datos:

> kernel.sf = st_as_sf(mydata.kernel.poly)
> st_intersection(kernel.sf, pa_data)
Simple feature collection with 0 features and 28 fields

devuelve 0 características porque hay cero intersección. Si hubiera una intersección, se podría utilizar st_area para obtener el área de intersección y calcular el porcentaje del polígono del hábitat en el área protegida.