In albhasan/sits.prodes: What the Package Does (One Line, Title Case)

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

suppressMessages(library(dplyr))
suppressMessages(library(ggplot2))
suppressMessages(library(kableExtra))
suppressMessages(library(knitr))
suppressMessages(library(magrittr))
suppressMessages(library(purrr))
suppressMessages(library(tidyr))
suppressMessages(library(sf))
suppressMessages(library(sits.prodes))
suppressMessages(library(xtable))

utils::data("expert_validation", package = "sits.prodes")

# Get the accuracies out of a confusionMatrix object. 
# 
# @param  x A confusionMatrix object (see caret::confusionMatrix).
# @return A list of numeric. The overall, user, and producer accuracies.
.get_accuracies <- function(x){
    if (length(x) == 1 && is.na(x)) 
        return(tibble::tibble(NA))
    res <- sits.prodes::asses_accuracy_simple(x$table)
    c(res$user, res$producer) %>%
        tibble::enframe(name = NULL) %>%
        dplyr::bind_cols(cname =  c(paste0("ua_", names(res$user)),
                                    paste0("pa_", names(res$producer)))) %>%
        tidyr::spread("cname", "value") %>%
        return()
}


# Check if the confusion matrix's values are greater than zero. 
# 
# @param  x         A confusionMatrix object (see caret::confusionMatrix).
# @param  label_vec A chatacter. Names of the labels to consider.
# @return           A length-one logical.
.is_mt_complete <- function(x, label_vec){
    if (length(x) == 1 && is.na(x)) 
        return(FALSE)
    mt <- x$table %>% as.matrix()
    if (!all(label_vec %in% colnames(mt)) && all(label_vec %in% rownames(mt)))
        return(FALSE)
    mt_lab <- mt[label_vec, label_vec]
    if (any(is.na(mt_lab)))
        return(FALSE)
    if (any(apply(mt_lab, 2, sum) == 0))
        return(FALSE)
    return(TRUE)
}

PRODES mapping to English

prodes_labels %>% 
    dplyr::select(-id_pd) %>% 
    knitr::kable() %>% 
    kableExtra::kable_styling(full_width = F)

Results

PRODES is the reference map (it is represented in the columns of the confusion matrices).
The data given in numer of pixels (except for user and producer accuracies).
Experts validation is a set of sample points for validating classification results. The samples were selected randomically an classified by hand by an expert (Rodrigo Begotti) using Landsat mosaics from PRODES Digital.

# get data for plots
acc_tb <- expert_validation %>%
    dplyr::mutate(
        is_mt_complete = purrr::map_lgl(.$confusion_matrix, .is_mt_complete, 
                                    label_vec = c("deforestation", "forest")),
        up_acc2 = purrr::map(.$confusion_matrix, .get_accuracies)
    ) %>%
    dplyr::filter(is_mt_complete == TRUE) %>%
    tidyr::drop_na(confusion_matrix) %>% 
    dplyr::select(experiment, algorithm, smooth, scene, pyear, up_acc2) %>%
    tidyr::unnest() %>%
    dplyr::select(-tidyselect::ends_with(match = "no forest"), 
                  -tidyselect::ends_with(match = "water")) %>%
    tidyr::gather(tidyselect::starts_with("pa_"), tidyselect::starts_with("ua_"),
                  key = "variable", value = "value") %>%
    dplyr::mutate(acc_type = ifelse(stringr::str_sub(variable, 1, 3) == "ua_",
                                    "user", "producer"),
                  variable = stringr::str_sub(variable, 4))

acc_pr <- acc_tb %>% dplyr::filter(acc_type == "producer") %>%
    dplyr::rename(prod_acc = value) %>%
    dplyr::select(-acc_type)
acc_ur <- acc_tb %>% dplyr::filter(acc_type == "user") %>%
    dplyr::rename(user_acc = value) %>%
    dplyr::select(-acc_type)
pa_acc_tb <- acc_pr %>%
    dplyr::inner_join(acc_ur, by = c("experiment", "algorithm", "smooth", 
                                     "scene", "pyear", "variable"))

f_plot <- function(data_tb, my_algorithm, my_smooth){
  prodes_years <- data_tb %>% 
    dplyr::pull(pyear) %>% 
    unique() %>% 
    sort() 
  data_tb <- data_tb %>%
    dplyr::filter(smooth == !!my_smooth, 
                  algorithm == !!my_algorithm)
  if (nrow(data_tb) < 1 || ncol(data_tb) < 1) {
    warning(sprintf("No data found for  %s %s", my_algorithm, my_smooth))
    return()
  }
  plot(
    data_tb %>% 
      ggplot2::ggplot(aes(x = prod_acc, y = user_acc,
                          shape = variable, color = experiment)) +
      ggplot2::geom_point(size = 2) +
      ggplot2::coord_fixed() +
      ggplot2::xlim(0, 1) +
      ggplot2::ylim(0, 1) +
      ggplot2::guides(shape = guide_legend(override.aes = list(size = 3))) +
      ggplot2::facet_wrap(scene ~ pyear, ncol = length(prodes_years)) +
      ggplot2::labs(title = sprintf("Expert validation %s %s", 
                                    my_algorithm, my_smooth),
                    x = "Producer accuracy",
                    y = "User accuracy")
  )
}

for (my_algorithm in unique(pa_acc_tb$algorithm)) {
    cat("\n\n")
    cat(paste0("#### ", my_algorithm), "\n")
    for (my_smooth in unique(pa_acc_tb$smooth)) {
        print(paste(my_algorithm, my_smooth, sep = "-"))
      #-----
      #### results_vote 
#[1] "results_vote-smooth_3x3_n10"
#[1] "results_vote-smooth_5x5_n10"
#[1] "results_vote-smooth_7x7_n10"
#[1] "results_vote-smooth_9x9_n10"
#[1] "results_vote-smooth_11x11_n10"
#[1] "results_vote-"
        suppressWarnings(f_plot(data_tb = pa_acc_tb, my_algorithm, my_smooth))
    }
}

for (i in 1:nrow(expert_validation)) {
    cap <- sprintf("Confusion matrix %s %s", unlist(expert_validation[i, "scene"]), 
                   unlist(expert_validation[i, "pyear"]))
    cat(paste0("### ", cap), "\n")
    if (!(length(expert_validation$confusion_matrix[[i]]) == 1 &&
          is.na(expert_validation$confusion_matrix[[i]]))) {
        p_conmat <- expert_validation$confusion_matrix[[i]]$table %>%
            apply(2, function(x){x/1}) %>% 
            add_upacc() %>%
            tibble::as_tibble(rownames = NA, caption = cap) %>% 
            knitr::kable(digits = 2, row.names = TRUE) %>% 
            kableExtra::kable_styling(full_width = F)
    }
    print(p_conmat)
    cat("\n\n")
}