Philippines LULC data

Code

library(terra)
library(sf)
library(ggplot2)
library(dplyr)
library(data.table)
library(knitr)

Code

area <- list.files("data/umrbpl/", full.names = TRUE, pattern = "shp$") |>
  read_sf() |>
  st_transform(crs = "EPSG:4326")
tmf2020 <- "data/tmf/JRC_TMF_AnnualChange_v1_2020_ASI_ID76_N20_E120.tif" |>
  rast() |>
  crop(area) |>
  mask(area)

Code

names_class <- c("LC1988", "LC2003", "AGG14", "AGG12", "CLASS_NAME")
lulc <- list.files("data/LULC", recursive = TRUE, "shp", full.names = TRUE) |>
  grep(pattern = "UMRBPL", value=TRUE) |>
  lapply(vect) |>
  # project in longlat
  lapply(function(x) project(x, "epsg:4326")) |>
  # harmonize class names
  lapply(function(x) {
    x$LC <- as.data.frame(x)[, which(names(x) %in% names_class)[1]]
    x$LC[grepl("wood|shrub", tolower(x$LC))] <- "Shrubs"
    x$LC[grepl("crop|grass|built|water|barren|fish", tolower(x$LC))] <- "Other"
    return(x) # nolint
  }) |>
  lapply(function(x) {
    rasterize(x, tmf2020, "LC")
  }) |>
  rast()

Temporal trends

Code

df_lulc <- as.data.frame(lulc, xy = TRUE)
names(df_lulc) <-
  c("lon", "lat", "LC1988", "LC2003", "AGG14", "AGG12", "CLASS_NAME")
df_lulc |>
  setDT() |>
  melt(id.vars = c("lon", "lat")) |>
  ggplot(aes(lon, lat, fill = value)) +
  geom_raster() +
  scale_fill_manual(values = c("chartreuse4", "chartreuse2", "grey", "gold")) +
  facet_wrap(~variable) +
  coord_equal() +
  theme_classic() +
  labs(x = NULL, y = NULL, fill = NULL) +
  theme(legend.position = c(0.85, 0.2))

Code

files <- list.files("data/LULC", recursive = TRUE, "shp") |>
  grep(pattern = "UMRBPL", value=TRUE)
years <- data.frame(
  layer = seq(1, length(files)),
  year = as.numeric(substr(files, 1, 4))
)
freq(lulc) |>
  mutate(area = count * 0.09) |>
  left_join(years) |>
  ggplot(aes(year, area, col = value)) +
  geom_line() +
  scale_color_manual(values = c("chartreuse4", "chartreuse2", "grey", "gold")) +
  theme_classic() +
  labs(x = "Year", y = "Area [ha]", col = NULL)

Comparison with TMF

Code

levels(tmf2020) <- data.frame(
  value = 0:6,
  category = c(
    NA, "Undisturbed Forest", "Degraded Forest", "Deforested",
    "Regrowing Forest", "Water", "Other"
  )
)

as.data.frame(c(tmf2020, lulc[[5]])) |>
  select(LC, category) |>
  table() |>
  kable()

	Undisturbed Forest	Degraded Forest	Deforested	Regrowing Forest	Water	Other
Closed Forest	25515	253	4	23	0	3
Open Forest	39846	5957	873	2188	0	829
Other	668	2240	4453	12158	33	55434
Shrubs	5275	15854	12885	47470	1	69921

Code

as.data.frame(c(tmf2020, lulc[[5]]), xy = TRUE) |>
  filter(!is.na(LC)) |>
  ggplot(aes(x, y, fill = category)) +
  geom_raster() +
  coord_equal() +
  theme_classic() +
  facet_wrap(~LC) +
  labs(x = NULL, y = NULL, fill = NULL)