R/sits_summary.R

In e-sensing/sits: Satellite Image Time Series Analysis for Earth Observation Data Cubes

#' @title  Summarize sits
#' @method summary sits
#' @name summary.sits
#' @author Gilberto Camara, \email{gilberto.camara@@inpe.br}
#' @author Felipe Souza, \email{felipe.souza@@inpe.br}
#' @description This is a generic function. Parameters depend on the specific
#' type of input.
#'
#' @param  object      Object of classes "sits".
#' @param  ...         Further specifications for \link{summary}.
#'
#' @return A summary of the sits tibble.
#'
#' @examples
#' if (sits_run_examples()) {
#'     summary(samples_modis_ndvi)
#' }
#'
#' @export
summary.sits <- function(object, ...) {
    # get frequency table
    data_labels <- table(object$label)

    # compose tibble containing labels, count and relative frequency columns
    result <- tibble::as_tibble(list(
        label = names(data_labels),
        count = as.integer(data_labels),
        prop = as.numeric(prop.table(data_labels))
    ))
    return(result)
}

#' @title  Summarize accuracy matrix for training data
#' @method summary sits_accuracy
#' @name summary.sits_accuracy
#' @author Gilberto Camara, \email{gilberto.camara@@inpe.br}
#' @description This is a generic function. Parameters depend on the specific
#' type of input.
#'
#' @param  object      Object of classe "sits_accuracy".
#' @param  ...         Further specifications for \link{summary}.
#'
#' @return A summary of the sample accuracy
#'
#' @examples
#' if (sits_run_examples()) {
#'     data(cerrado_2classes)
#'     # split training and test data
#'     train_data <- sits_sample(cerrado_2classes, frac = 0.5)
#'     test_data  <- sits_sample(cerrado_2classes, frac = 0.5)
#'     # train a random forest model
#'     rfor_model <- sits_train(train_data, sits_rfor())
#'     # classify test data
#'     points_class <- sits_classify(
#'         data = test_data,
#'         ml_model = rfor_model
#'     )
#'     # measure accuracy
#'     acc <- sits_accuracy(points_class)
#'     summary(acc)
#' }
#'
#' @export
summary.sits_accuracy <- function(object, ...) {
    sits_accuracy_summary(object)
}
#' @title  Summarize accuracy matrix for area data
#' @method summary sits_area_accuracy
#' @name summary.sits_area_accuracy
#' @author Gilberto Camara, \email{gilberto.camara@@inpe.br}
#' @description This is a generic function. Parameters depend on the specific
#' type of input.
#'
#' @param  object      Object of classe "sits_accuracy".
#' @param  ...         Further specifications for \link{summary}.
#'
#' @return A summary of the sample accuracy
#'
#' @examples
#' if (sits_run_examples()) {
#'     # create a data cube from local files
#'     data_dir <- system.file("extdata/raster/mod13q1", package = "sits")
#'     cube <- sits_cube(
#'         source = "BDC",
#'         collection = "MOD13Q1-6",
#'         data_dir = data_dir
#'     )
#'     # create a random forest model
#'     rfor_model <- sits_train(samples_modis_ndvi, sits_rfor())
#'     # classify a data cube
#'     probs_cube <- sits_classify(
#'         data = cube, ml_model = rfor_model, output_dir = tempdir()
#'     )
#'     # label the probability cube
#'     label_cube <- sits_label_classification(
#'         probs_cube,
#'         output_dir = tempdir()
#'     )
#'     # obtain the ground truth for accuracy assessment
#'     ground_truth <- system.file("extdata/samples/samples_sinop_crop.csv",
#'         package = "sits"
#'     )
#'     # make accuracy assessment
#'     as <- sits_accuracy(label_cube, validation = ground_truth)
#'     summary(as)
#' }
#'
#' @export
summary.sits_area_accuracy <- function(object, ...) {
    print.sits_area_accuracy(object)
}
#' @title  Summarize data cubes
#' @method summary raster_cube
#' @name summary.raster_cube
#' @author Gilberto Camara, \email{gilberto.camara@@inpe.br}
#' @author Felipe Souza, \email{felipe.souza@@inpe.br}
#' @description This is a generic function. Parameters depend on the specific
#' type of input.
#'
#' @param  object      Object of classes "raster_cube".
#' @param  ...         Further specifications for \link{summary}.
#' @param  tile        Tile to be summarized
#' @param  date        Date to be summarized
#'
#' @return A summary of the data cube.
#'
#' @examples
#' if (sits_run_examples()) {
#'     # create a data cube from local files
#'     data_dir <- system.file("extdata/raster/mod13q1", package = "sits")
#'     cube <- sits_cube(
#'         source = "BDC",
#'         collection = "MOD13Q1-6",
#'         data_dir = data_dir
#'     )
#'     summary(cube)
#' }
#'
#' @export
summary.raster_cube <- function(object, ..., tile = NULL, date = NULL) {
    # Pre-conditional check
    .check_date_parameter(date, allow_null = TRUE)
    .check_chr_parameter(tile, allow_null = TRUE)
    # Extract the chosen tile
    if (!is.null(tile)) {
        object <- .summary_check_tile(object, tile)
    }
    # Extract the chosen date
    if (!is.null(date)) {
        object <- .cube_filter_dates(object, dates = date)
    }
    # Display cube general metadata
    cli::cli_h1("Cube Metadata")
    cli::cli_li("Class: {.field raster_cube}")
    cube_bbox <- sits_bbox(object)[, c('xmin', 'xmax', 'ymin', 'ymax')]
    cli::cli_li("Bounding Box: xmin = {.field {cube_bbox[['xmin']]}},
                               xmax = {.field {cube_bbox[['xmax']]}},
                               ymin = {.field {cube_bbox[['ymin']]}},
                               ymax = {.field {cube_bbox[['ymax']]}}")
    cli::cli_li("Bands: {.field {sits_bands(object)}}")
    timeline <- unique(lubridate::as_date(unlist(.cube_timeline(object))))
    cli::cli_li("Timeline: {.field {timeline}}")
    is_regular <- .cube_is_regular(object)
    cli::cli_li("Regular cube: {.field {is_regular}}")
    # Display cube cloud coverage
    if ("CLOUD" %in% .cube_bands(object) &&
        .has_column(.fi(object), "cloud_cover")) {
        cube_unnest <- tidyr::unnest(
            object[, c("tile", "file_info")], "file_info"
        )
        cli::cli_h1("Cloud cover info")
        cube_unnest <- cube_unnest[, c("tile", "date", "cloud_cover")]
        cube_unnest <- unique(dplyr::arrange(cube_unnest, .data[["date"]]))
        print(cube_unnest, n = Inf)
    }
    # Display raster summary
    cli::cli_h1("Cube Summary")
    sum <- slider::slide(object, function(tile) {
        # Get the first date to not read all images
        date <- .default(date, .tile_timeline(tile)[[1]])
        tile <- .tile_filter_dates(tile, date)
        bands <- if (is_regular) .tile_bands(tile) else .tile_bands(tile)[[1]]
        tile <- .tile_filter_bands(tile, bands)
        cli::cli_h3("Tile: {.field {tile$tile}} and Date: {.field {date}}")
        rast <- .raster_open_rast(.tile_paths(tile))
        sum <- suppressWarnings(.raster_summary(rast))
        print(sum)
        return(sum)
    })
    # Return the summary from the cube
    names(sum) <- .cube_tiles(object)
    return(invisible(sum))
}
#' @title Summary of a derived cube
#' @author Felipe Souza, \email{felipe.souza@@inpe.br}
#' @noRd
#' @param object data cube
#' @param  ...         Further specifications for \link{summary}.
#' @param tile A \code{tile}.
#' @return Summary of a derived cube
#'
#' @examples
#' if (sits_run_examples()) {
#'     # create a data cube from local files
#'     data_dir <- system.file("extdata/raster/mod13q1", package = "sits")
#'     cube <- sits_cube(
#'         source = "BDC",
#'         collection = "MOD13Q1-6",
#'         data_dir = data_dir
#'     )
#'     # create a random forest model
#'     rfor_model <- sits_train(samples_modis_ndvi, sits_rfor())
#'     # classify a data cube
#'     probs_cube <- sits_classify(
#'         data = cube, ml_model = rfor_model, output_dir = tempdir()
#'     )
#'     summary(probs_cube)
#'     # get the variance cube
#'     variance_cube <- sits_variance(
#'         probs_cube,
#'         output_dir = tempdir()
#'     )
#'     summary(variance_cube)
#' }
#'
#' @export
summary.derived_cube <- function(object, ..., tile = NULL) {
    # Pre-conditional check
    .check_chr_parameter(tile, allow_null = TRUE)
    # Extract the chosen tile
    if (!is.null(tile)) {
        object <- .summary_check_tile(object, tile)
    }
    # Display cube general metadata
    cli::cli_h1("Cube Metadata")
    cli::cli_li("Class: {.field derived_cube}")
    cube_bbox <- sits_bbox(object)[, c('xmin', 'xmax', 'ymin', 'ymax')]
    cli::cli_li("Bounding Box: xmin = {.field {cube_bbox[['xmin']]}},
                               xmax = {.field {cube_bbox[['xmax']]}},
                               ymin = {.field {cube_bbox[['ymin']]}},
                               ymax = {.field {cube_bbox[['ymax']]}}")
    cli::cli_li("Band(s): {.field {sits_bands(object)}}")
    timeline <- unique(lubridate::as_date(unlist(.cube_timeline(object))))
    cli::cli_li("Timeline: {.field {timeline}}")
    # get sample size
    sample_size <- .conf("summary_sample_size")
    # Get tile name
    tile <- .default(tile, .cube_tiles(object)[[1]])
    cli::cli_h1("Cube Summary")
    tile <- .cube_filter_tiles(object, tile)
    # get the bands
    band <- sits_bands(tile)
    .check_num(
        x = length(band),
        min = 1,
        max = 1,
        is_integer = TRUE,
        msg = "invalid cube - more than one probs band"
    )
    # extract the file paths
    files <- .tile_paths(tile)
    # read the files with terra
    r <- terra::rast(files)
    # get the a sample of the values
    values <- r |>
        terra::spatSample(size = sample_size, na.rm = TRUE)
    # scale the values
    band_conf <- .tile_band_conf(tile, band)
    scale <- .scale(band_conf)
    offset <- .offset(band_conf)
    sum <- summary(values * scale + offset)
    colnames(sum) <- sits_labels(tile)
    return(sum)
}
#' @title  Summarize data cubes
#' @method summary class_cube
#' @name summary.class_cube
#' @author Gilberto Camara, \email{gilberto.camara@@inpe.br}
#' @description This is a generic function. Parameters depend on the specific
#' type of input.
#' @param  object    Object of class "class_cube"
#' @param ...        Further specifications for \link{summary}.
#' @param  tile        Tile to be summarized
#'
#' @return A summary of a classified cube
#'
#' @examples
#' if (sits_run_examples()) {
#'     # create a data cube from local files
#'     data_dir <- system.file("extdata/raster/mod13q1", package = "sits")
#'     cube <- sits_cube(
#'         source = "BDC",
#'         collection = "MOD13Q1-6",
#'         data_dir = data_dir
#'     )
#'     # create a random forest model
#'     rfor_model <- sits_train(samples_modis_ndvi, sits_rfor())
#'     # classify a data cube
#'     probs_cube <- sits_classify(
#'         data = cube, ml_model = rfor_model, output_dir = tempdir()
#'     )
#'     # label the probability cube
#'     label_cube <- sits_label_classification(
#'         probs_cube,
#'         output_dir = tempdir()
#'     )
#'     summary(label_cube)
#' }
#' @export
#'
summary.class_cube <- function(object, ..., tile = NULL) {
    # Pre-conditional check
    .check_chr_parameter(tile, allow_null = TRUE)
    # Extract the chosen tile
    if (!is.null(tile)) {
        object <- .summary_check_tile(object, tile)
    }
    # Display cube general metadata
    cli::cli_h1("Cube Metadata")
    cli::cli_li("Class: {.field class_cube}")
    cube_bbox <- sits_bbox(object)[, c('xmin', 'xmax', 'ymin', 'ymax')]
    cli::cli_li("Bounding Box: xmin = {.field {cube_bbox[['xmin']]}},
                               xmax = {.field {cube_bbox[['xmax']]}},
                               ymin = {.field {cube_bbox[['ymin']]}},
                               ymax = {.field {cube_bbox[['ymax']]}}")
    cli::cli_li("Band(s): {.field {sits_bands(object)}}")
    timeline <- unique(lubridate::as_date(unlist(.cube_timeline(object))))
    cli::cli_li("Timeline: {.field {timeline}}")
    # Get tile name
    tile <- .default(tile, .cube_tiles(object)[[1]])
    cli::cli_h1("Cube Summary")
    tile <- .cube_filter_tiles(object, tile)
    # get the bands
    band <- sits_bands(tile)
    .check_num(
        x = length(band),
        min = 1,
        max = 1,
        is_integer = TRUE,
        msg = "invalid cube - more than one probs band"
    )
    # extract the file paths
    files <- .tile_paths(tile)
    # read raster files
    r <- .raster_open_rast(files)
    # get a frequency of values
    class_areas <- .raster_freq(r)
    # transform to km^2
    cell_size <- .tile_xres(tile) * .tile_yres(tile)
    class_areas[["area"]] <-  (class_areas[["count"]] * cell_size) / 10^6
    # change value to character
    class_areas <- dplyr::mutate(class_areas,
        value = as.character(.data[["value"]])
    )
    # create a data.frame with the labels
    labels <- sits_labels(tile)
    df1 <- data.frame(value = names(labels), class = unname(labels))
    # join the labels with the areas
    sum <- dplyr::full_join(df1, class_areas, by = "value")
    sum <- dplyr::mutate(sum,
        area_km2 = signif(.data[["area"]], 3),
        .keep = "unused"
    )
    # remove layer information
    sum_clean <- sum[, -3] |>
        stats::na.omit()
    # are there NA's ?
    sum_NA <- dplyr::filter(sum, is.na(.data[["area_km2"]]))
    # inform missing classes
    if (nrow(sum_NA) > 0) {
        message(
            "classes ", paste0(sum_NA$class, collapse = " "),
            " have no area"
        )
    }
    # show the result
    return(sum_clean)
}