Land Cover Change with rgee

Land Cover change with Corine CLC

# packages
pacman::p_load('ps', "rgee", "tidyverse", "sf", "ggplot2", "patchwork",
               'reticulate', 'googledrive', 'stars', 'plotKML', 'readxl',
               'networkD3', 'OpenLand', 'ggtern')

ee_Initialize(email = '[email protected]')

# CORINE
clc18 = ee$Image('COPERNICUS/CORINE/V20/100m/2018')$select('landcover');
clc12 = ee$Image('COPERNICUS/CORINE/V20/100m/2012')$select('landcover');
clc1812 = clc18$multiply(1000)$add(clc12);

# CLC colour scheme
# read corine color palette
colours <- readxl::read_xls('D:/Dropbox/programacao/gee/clc2000legend.xls')
colourshex <- as_tibble(str_split_fixed(colours$RGB, "-", 3)) %>% mutate_if(is.character, as.integer)
colourshex$CLC_CODE <- colours$CLC_CODE
colourshex <- colourshex %>%  filter(complete.cases(.)) %>%  mutate(hexcode = rgb2hex(V1, V2, V3))

# ROI
roi <- st_buffer(sf::st_sfc(st_point(c(-8.9373, 38.8900))), 0.1) %>% sf::st_set_crs(4326)
roi_ee <- roi %>% sf_as_ee()

# Map View
Map$setCenter(-8.94, 38.89, 10)
Map$addLayer(clc18, {}, 'Land Cover 18', opacity = 0.4) +
  Map$addLayer(roi_ee, {}, 'roi',  opacity = 0.4)

# Reduce to region
vListchange <- tibble(change = clc1812$reduceRegion(
  reducer= ee$Reducer$toList(),
  maxPixels = 1e9,
  geometry =  roi_ee
)$values()$get(0)$getInfo()
)

# Summarise
change1812 <- vListchange %>% separate(change, into = c('clc18', 'clc12'), sep = 3, remove = F)
change1812un <- distinct(change1812)
legend_tab <- tibble(categoryValue = as.numeric(unique(c(change1812un$clc18,change1812un$clc18)))) %>%
  left_join(select(colourshex, CLC_CODE, color= hexcode), by = c('categoryValue' = 'CLC_CODE')) %>%
  mutate(categoryName=as.factor(categoryValue)) %>%
  select(categoryValue, categoryName, color)

summchange1812 <- change1812 %>% group_by(change) %>%
  summarise(area_ha = n(), .groups = 'drop') %>%
  inner_join(distinct(change1812), by = c('change' = 'change')) %>%
  mutate(changed = clc18 != clc12,
         CLC_CODE = as.numeric(clc18),
         fromto = paste(clc12, clc18, sep = '-')) %>%
  inner_join(select(colourshex, CLC_CODE, hexcode), by = c('CLC_CODE' = 'CLC_CODE'))

# ChordDiagram inputs
conting_tab <- tibble(Period = '2012-2018',
                      From = as.integer(summchange1812$clc12),
                      To = as.integer(summchange1812$clc18),
                      km2 = summchange1812$area_ha/1000,
                      QtPixel = summchange1812$area_ha,
                      Interval = 6,
                      yearFrom = 2012,
                      yearTo = 2018)

p1 <- chordDiagramLand(dataset = conting_tab, legendtable = legend_tab, legtitle = "CLC Class",
                       legendsize = .8)