<jupyter_start><jupyter_text>Land use and land cover classification in the Brazilian Cerrado biome using Brazil Data Cube    <NAME>  , <NAME>  , <NAME>  , <NAME>  ,    <NAME>  , <NAME>  , <NAME>*          Earth Observation and Geoinformatics Division, National Institute for Space Research (INPE)        Avenida dos Astronautas, 1758, Jardim da Granja, São José dos Campos, SP 12227-010, Brazil        * Author to whom correspondence should be addressed.        February 24, 2021Abstract. This Jupyter Notebook compendium contains useful information for the creation of land use and land cover (LULC) maps using Earth observations data cubes and machine learning (ML) techniques. The code is based on the research pipeline described in the paper Earth Observation Data Cubes for Brazil: Requirements, Methodology and Products. These notebooks access open data available in the Brazil Data Cube platform.        This Jupyter Notebook is supplement to the Section 5 of the following paper:        <NAME>.; <NAME>.; <NAME>.; <NAME>.; <NAME>.; <NAME>.; <NAME>.; <NAME>.; <NAME>. 2020. Earth Observation Data Cubes for Brazil: Requirements, Methodology and Products. Remote Sens. 12, no. 24: 4033. DOI: 10.3390/rs12244033.      Land Use and Cover Mapping from CBERS-4/AWFI Data Cubes This document will present the steps used to generate the CBERS-4/AWFI classification map presented in the paper. As presented in the article, the classification process was done using the [SITS R package](https://github.com/e-sensing/sits).  Study Area and samples In this application we use the same region of interest and samples described in [Ferreira et al. (2020)](https://doi.org/10.3390/rs12244033). As depicted in Figure 1, the region is located in the Bahia state (Brazil), between the Cerrado and Caatinga biomes.  Figure 1 - Study area in relation to Brazil and its biomes.  Set a pseudo-randomic seed We will fix a pseudo-randomic seed in order to run the code:<jupyter_code>set.seed(777)<jupyter_output><empty_output><jupyter_text>Loading the software packages<jupyter_code>library(sits)<jupyter_output><empty_output><jupyter_text>The user should also provides his access key to the Brazil Data Cube platform:<jupyter_code>MY_ACCESS_KEY <- "change-me"
Sys.setenv(BDC_ACCESS_KEY = MY_ACCESS_KEY)<jupyter_output><empty_output><jupyter_text>Defining the Data Cube Let's start by defining the region of interest (ROI) as a sub-space of the red rectangle shown in Figure 1. The ROI is available in a file name `roi.rds` under the directory `roi`:<jupyter_code>roi <- readRDS(url("https://brazildatacube.dpi.inpe.br/geo-knowledge-hub/bdc-article/roi/roi.rds"))<jupyter_output><empty_output><jupyter_text>> The `roi` is a list with two components:> * `classification_roi`: contains the geometry boundary for the classification.> * `search_roi`: a smaller rectangle than the `classification_roi`, that intersects only the data cube tiles we are interested to use in the classification.Next we define a time interval based on the crop calendar to define the working period:<jupyter_code>start_date  <- "2018-09-01"
end_date    <- "2019-08-31"<jupyter_output><empty_output><jupyter_text>In this Jupyter Notebook we focus the classification based on a CBERS-4/AWFI data cube named `CB4_64_16D_STK-1`:<jupyter_code>collection  <- "CB4_64_16D_STK-1"<jupyter_output><empty_output><jupyter_text>Finally, let's define the data cube.The `sits` package will access the CBERS-4/AWFI data cube available in the Brazil Data Cube platform through the STAC web service:<jupyter_code>cube <- sits_cube(
  type        = "BDC",
  name        = "cube_to_classify",
  url         = "https://brazildatacube.dpi.inpe.br/stac/",
  collection  = collection,
  start_date  = start_date,
  end_date    = end_date,
  roi         = roi$search_roi
)<jupyter_output><empty_output><jupyter_text>> The definition above includes the spectral bands `Red`, `Green`, `Blue`, `Near-Infrared (NIR)` and the vegetation indices `EVI` and `NDVI` already available in the cube.> It also limits the temporal extension to `2018-09` to `2019-08`.  Loading the Training Samples Now, let's load the samples from a prepared file named `CB4_64_16D_STK_1.rds`:<jupyter_code>samples <- readRDS(url("https://brazildatacube.dpi.inpe.br/geo-knowledge-hub/bdc-article/training-samples/rds/CB4_64_16D_STK_1.rds"))<jupyter_output><empty_output><jupyter_text>> The Jupyter Notebook entitled [Extracting time series from sample locations](./01_ExtractTimeSeries.ipynb) describes in detail how to prepare this file.  MultiLayer Perceptron model definition For the classification of data cubes, the article presents the use of an MLP network with five hidden layers with 512 neurons, trained with the backpropagation algorithm, using the Adam optimizer. The model uses the ReLu activation function.Below is the definition of this model using the [SITS package](https://github.com/e-sensing/sits).<jupyter_code>mlp_model <- sits_deeplearning(layers        = c(512, 512, 512, 512, 512),
                               activation    = "relu",
                               optimizer     = keras::optimizer_adam(lr = 0.001),
                               epochs        = 200)<jupyter_output><empty_output><jupyter_text>Below, the defined model is trained using the same samples used in the article.<jupyter_code>dl_model <- sits_train(samples, mlp_model)<jupyter_output><empty_output><jupyter_text>Output Directory All the results generated in this document will be saved in your user's home directory, inside `bdc-article/results/CB4_64_16D_STK_1` directory:<jupyter_code>output_dir <- "results/CB4_64_16D_STK_1"

dir.create(
    path         = output_dir, 
    showWarnings = FALSE, 
    recursive    = TRUE
)<jupyter_output><empty_output><jupyter_text>Classifying the tiles from the data cube Before running the classification step, you should define the hardware resources that `sits` will be allowed to use during the classification:<jupyter_code>classification_memsize    <- 8 # GB
classification_multicores <- 3 # CPU logical cores<jupyter_output><empty_output><jupyter_text>The next cell uses the trainned MLP model (`dl_model`) to perform a classification based on the temporal data from the data cube: > This is a time-consuming process.<jupyter_code>probs <- sits_classify(data       = cube,
                       ml_model   = dl_model,
                       memsize    = classification_memsize,
                       multicores = classification_multicores,
                       roi        = roi$classification_roi,
                       output_dir = output_dir)<jupyter_output><empty_output><jupyter_text>> The generated data will be stored under the directory indicated by `output_dir`.> Note that here we use a geometry boundary from `roi$classification_roi` that is smaller than the region defined by the samples.The classification output of the Multilayer Perceptron (MLP) model is a raster with three layers (one for each land use class) containing the probabilities of each pixel belonging to each of the classes.The raster file named `cube_to_classify_022024_probs_2018_8_2019_7_v1.tif` has 3 layers containing scaled probabilities (`x 10,000`) corresponding to the classes `Crop` (layer 1), `Natural Vegetation` (layer 2), and `Pasture` (layer 3).  Generating the Thematic Map We are going to apply a probability Bayesian smoother method over the output of the MLP. This procedure uses the information of a pixel’s neighborhood to update its probabilities by taking the maximum likelihood estimator. The smoothing procedure removes isolated pixel class values and produces more homogeneous spatial areas.The next cell perform this operation:<jupyter_code>probs_smoothed <- sits_smooth(probs, type = "bayes", output_dir = output_dir)<jupyter_output><empty_output><jupyter_text>> The above cell will create a file name `cube_to_classify_022024_probs_2018_8_2019_7_bayes_v1.tif`.After that, to generate the thematic maps the most probable class is taken as the pixel class. Each class is represented by the codes 1 (Crop), 2 (Natural Vegetation), and 3 (Pasture). The next cell show how to perform this step:<jupyter_code>labels <- sits_label_classification(probs_smoothed, output_dir = output_dir)<jupyter_output><empty_output><jupyter_text>> The final map is named `cube_to_classify_022024_probs_class_2018_8_2019_7_v1.tif`.  Visualizing the Thematic Map Finally, let's use the `rgdal` library to plot the resulted map:<jupyter_code>library(rgdal)
rst <- raster::raster(
    paste0(output_dir, "/cube_to_classify_022024_probs_class_2018_8_2019_7_v1.tif")
)
plot(rst)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.07740367927802846, 'fraction_numerical': 0.028578040032396158, 'mean_word_length': 4.6202860858257475, 'pattern_counts': {'":': 0, '<': 97, '<?xml version=': 0, '>': 92, 'https://': 6, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '16', 'dir': 'jupyter-structured-clean-dedup', 'id': '71468', 'license': 'MIT', 'n_tokens_mistral': 2797, 'n_tokens_neox': 2534, 'n_words': 919, 'path': 'jupyter/R/bdc-article/02_CB4_64_16D_STK-1_Classification.ipynb', 'repo_name': 'lorenalves/code-gallery', 'size': '8940'}