Deer data for KK - LS - TP area

Author

Virginia Morera-Pujol

SMARTDEER models

Predictions of the SMARTDEER models over the area of interest.

  • Values go from 0 to 1 nationwide, representing increasing likelihood of finding at least 1 deer in that area.

  • Resolution of the raster layer is 25 km2

  • Deer data used for the model span years 2007 - 2022

Reference: Morera-Pujol, V., Mostert, P.S., Murphy, K.J., Burkitt, T., Coad, B., McMahon, B.J., Nieuwenhuis, M., Morelle, K., Ward, A.I. and Ciuti, S. (2023), Bayesian species distribution models integrate presence-only and presence–absence data to predict deer distribution and relative abundance. Ecography, 2023: e06451. https://doi.org/10.1111/ecog.06451

Code 
fallow_IG <- project(fallow, crs("EPSG:29902"))
fallow_mean <- rescale0to1(fallow_IG$mean)

fallow_mask <- mask(crop(fallow_mean, counties), counties)

ggplot() +
  geom_spatraster(data = fallow_mask, aes(fill = mean)) +
  geom_sf(data = ded, fill = NA, col = "darkgray") +   
  geom_sf(data = ded2, col = "red", fill = NA) + 
  geom_sf(data = counties, col = "black", fill = NA) +
  scale_fill_viridis_c(na.value = NA) + 
  theme_bw() + 
  ggtitle("Fallow deer")

Code 
sika_IG <- project(sika, crs("EPSG:29902"))
sika_mean <- rescale0to1(sika_IG$mean)

sika_mask <- mask(crop(sika_mean, counties), counties)

ggplot() +
  geom_spatraster(data = sika_mask, aes(fill = mean)) +
  geom_sf(data = ded, fill = NA, col = "darkgray") + 
  geom_sf(data = ded2, col = "red", fill = NA) + 
  geom_sf(data = counties, col = "black", fill = NA) +
  scale_fill_viridis_c(na.value = NA) + 
  theme_bw() + 
  ggtitle("Sika deer")

Code 
red_IG <- project(red, crs("EPSG:29902"))
red_mean <- rescale0to1(red_IG$mean)

red_mask <- mask(crop(red_mean, counties), counties)

ggplot() +
  geom_spatraster(data = red_mask, aes(fill = mean)) +
  geom_sf(data = ded, fill = NA, col = "darkgray") + 
  geom_sf(data = ded2, col = "red", fill = NA) + 
  geom_sf(data = counties, col = "black", fill = NA) +
  scale_fill_viridis_c(na.value = NA) + 
  theme_bw() + 
  ggtitle("Red deer")

Disaggregation models

Predictions from the disaggregation models for 2018 and the area of interest

  • Resolution of the raster layer is 25 km2

  • Deer data used for the model span years 2007 - 2022

  • Data comes from the county-level culling returns and have been modelled to a higher resolution

  • Values represent the number of deer estimated to be culled in each 5 by 5 km cell

Reference: Murphy, Kilian J., Simone Ciuti, Tim Burkitt, and Virginia Morera-Pujol. 2023. “ Bayesian Areal Disaggregation Regression to Predict Wildlife Distribution and Relative Density with Low-Resolution Data.” Ecological Applications 33(8): e2924. https://doi.org/10.1002/eap.2924

Code 
fallow_IG <- project(fallow18, crs("EPSG:29902"))

fallow_mask <- mask(crop(fallow_IG, counties), counties)

ggplot() +
  geom_spatraster(data = fallow_mask, aes(fill = Fallow18_Prediction)) +
  geom_sf(data = ded, fill = NA, col = "darkgray") + 
  geom_sf(data = counties, col = "white", fill = NA) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  # scale_fill_viridis_c(na.value = NA, trans = "log") + 
  scale_fill_viridis_c(na.value = NA) + 
  theme_bw() + 
  ggtitle("Culled fallow deer in 2018")

Code 
sika_IG <- project(sika18, crs("EPSG:29902"))

sika_mask <- mask(crop(sika_IG, counties), counties)

ggplot() +
  geom_spatraster(data = sika_mask, aes(fill = Sika18_Prediction)) +
  geom_sf(data = ded, fill = NA, col = "darkgray") + 
  geom_sf(data = counties, col = "white", fill = NA) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  # scale_fill_viridis_c(na.value = NA, trans = "log") + 
  scale_fill_viridis_c(na.value = NA) + 
  theme_bw() + 
  ggtitle("Culled sika deer in 2018")

Code 
red_IG <- project(red18, crs("EPSG:29902"))

red_mask <- mask(crop(red_IG, counties), counties)

ggplot() +
  geom_spatraster(data = red_mask, aes(fill = Red18_Prediction )) +
  geom_sf(data = ded, fill = NA, col = "darkgray") + 
  geom_sf(data = counties, col = "white", fill = NA) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  # scale_fill_viridis_c(na.value = NA, trans = "log") + 
  scale_fill_viridis_c(na.value = NA) + 
  theme_bw() + 
  ggtitle("Culled red deer in 2018")

Using the disaggregation predictions, we can calculate the total number of deer culled in the area, as well as the minimum, mean, and maximum number of deer culled by 25 km2 cell.

Code 
counties <- counties %>% 
  mutate(
    sum_red = exact_extract(x = red_mask, counties, 'sum', progress = FALSE),
    mean_red = exact_extract(x = red_mask, counties, 'mean', progress = FALSE),
    min_red = exact_extract(x = red_mask, counties, 'min', progress = FALSE),
    max_red = exact_extract(x = red_mask, counties, 'max', progress = FALSE),
    sum_sika = exact_extract(x = sika_mask, counties, 'sum', progress = FALSE),
    mean_sika = exact_extract(x = sika_mask, counties, 'mean', progress = FALSE),
    min_sika = exact_extract(x = sika_mask, counties, 'min', progress = FALSE),
    max_sika = exact_extract(x = sika_mask, counties, 'max', progress = FALSE),
    sum_fallow = exact_extract(x = fallow_mask, counties, 'sum', progress = FALSE),
    mean_fallow = exact_extract(x = fallow_mask, counties, 'mean', progress = FALSE),
    min_fallow = exact_extract(x = fallow_mask, counties, 'min', progress = FALSE),
    max_fallow = exact_extract(x = fallow_mask, counties, 'max', progress = FALSE)
  )

# Fallow 
ggplot(counties) +
  geom_sf(aes(fill = sum_fallow)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Total culled animals", fill = "Fallow deer") +
  theme_bw() + 
  
  ggplot(counties) +
  geom_sf(aes(fill = mean_fallow)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Average per 25 sqkm", fill = "Fallow deer") +
  theme_bw() + 
  
  ggplot(counties) +
  geom_sf(aes(fill = min_fallow)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Minimum per 25 sqkm", fill = "Fallow of deer") +
  theme_bw() + 
  
  ggplot(counties) +
  geom_sf(aes(fill = max_fallow)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Maximum per 25 sqkm", fill = "Fallow deer") +
  theme_bw() + 
  
  plot_annotation(title = 'Fallow deer culling numbers', 
                  subtitle = "2894 deer culled in total in the area in 2018")

Code 
# Sika 
ggplot(counties) +
  geom_sf(aes(fill = sum_sika)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Total culled animals", fill = "sika deer") +
  theme_bw() + 
  
  ggplot(counties) +
  geom_sf(aes(fill = mean_sika)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Average per 25 sqkm", fill = "sika deer") +
  theme_bw() + 
  
  ggplot(counties) +
  geom_sf(aes(fill = min_sika)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Minimum per 25 sqkm", fill = "sika of deer") +
  theme_bw() + 
  
  ggplot(counties) +
  geom_sf(aes(fill = max_sika)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Maximum per 25 sqkm", fill = "sika deer") +
  theme_bw() + 
  
  plot_annotation(title = 'Sika deer culling numbers', 
                  subtitle = "507 deer culled in total in the area in 2018")

Code 
# Red 
ggplot(counties) +
  geom_sf(aes(fill = sum_red)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Total culled animals", fill = "Red deer") +
  theme_bw() + 
  
ggplot(counties) +
  geom_sf(aes(fill = mean_red)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Average per 25 sqkm", fill = "Red deer") +
  theme_bw() + 
  
ggplot(counties) +
  geom_sf(aes(fill = min_red)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Minimum per 25 sqkm", fill = "Red deer") +
  theme_bw() + 
  
ggplot(counties) +
  geom_sf(aes(fill = max_red)) +
  geom_sf(data = ded2, col = "red", fill = NA) + 
  scale_fill_viridis_c() + 
  labs(title = "Maximum per 25 sqkm", fill = "Red deer") +
  theme_bw() + 

plot_annotation(title = 'Red deer culling numbers', 
                subtitle = "429 deer culled in total in the area in 2018")