doc/tutorial.R

In floriandeboissieu/biodivMapR: biodivMapR: an R package for a- and ß-diversity mapping using remotely-sensed images

## ----setup, include = FALSE----------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  eval=FALSE
)

## ----Input / Output files------------------------------------------------
#  Input_Image_File  = system.file('extdata', 'RASTER', 'S2A_T33NUD_20180104_Subset', package = 'biodivMapR')
#  
#  # Input.Image.File  = raster2BIL(Raster.Path = Input.Image.File,
#  #                                        Sensor = 'SENTINEL_2A',
#  #                                        Convert.Integer = TRUE,
#  #                                        Output.Directory = '~/test')
#  
#  Input_Mask_File   = FALSE
#  
#  Output_Dir        = 'RESULTS'

## ----Spatial resolution--------------------------------------------------
#  window_size = 10

## ----PCA filtering-------------------------------------------------------
#  FilterPCA = FALSE

## ----Computing options---------------------------------------------------
#  nbCPU         = 2
#  MaxRAM        = 0.5
#  nbclusters    = 50

## ----Mask non vegetated / shaded / cloudy pixels-------------------------
#  NDVI_Thresh = 0.5
#  Blue_Thresh = 500
#  NIR_Thresh  = 1500
#  print("PERFORM RADIOMETRIC FILTERING")
#  Input_Mask_File = perform_radiometric_filtering(Input_Image_File, Input_Mask_File, Output_Dir,
#                                              NDVI_Thresh = NDVI_Thresh, Blue_Thresh = Blue_Thresh,
#                                              NIR_Thresh = NIR_Thresh)

## ----PCA-----------------------------------------------------------------
#  print("PERFORM PCA ON RASTER")
#  PCA_Output        = perform_PCA(Input_Image_File, Input_Mask_File, Output_Dir,
#                                 FilterPCA = TRUE, nbCPU = nbCPU,MaxRAM = MaxRAM)
#  # path for the PCA raster
#  PCA_Files         = PCA_Output$PCA_Files
#  # number of pixels used for each partition used for k-means clustering
#  Pix_Per_Partition = PCA_Output$Pix_Per_Partition
#  # number of partitions used for k-means clustering
#  nb_partitions     = PCA_Output$nb_partitions
#  # path for the updated mask
#  Input_Mask_File   = PCA_Output$MaskPath
#  # parameters of the PCA model
#  PCA_model         = PCA_Output$PCA_model
#  # definition of spectral bands to be excluded from the analysis
#  SpectralFilter    = PCA_Output$SpectralFilter
#  
#  print("Select PCA components for diversity estimations")
#  select_PCA_components(Input_Image_File, Output_Dir, PCA_Files, File_Open = TRUE)

## ----Spectral species map------------------------------------------------
#  print("MAP SPECTRAL SPECIES")
#  map_spectral_species(Input_Image_File, Output_Dir, PCA_Files, PCA_model, SpectralFilter, Input_Mask_File,
#                       Pix_Per_Partition, nb_partitions, nbCPU=nbCPU, MaxRAM=MaxRAM,
#                       nbclusters = nbclusters, TypePCA = TypePCA, CR = TRUE)

## ----alpha and beta diversity maps---------------------------------------
#  print("MAP ALPHA DIVERSITY")
#  # Index.Alpha   = c('Shannon','Simpson')
#  Index_Alpha   = c('Shannon')
#  map_alpha_div(Input_Image_File, Output_Dir, window_size, nbCPU=nbCPU,
#                MaxRAM=MaxRAM, Index_Alpha = Index_Alpha, nbclusters = nbclusters)
#  
#  print("MAP BETA DIVERSITY")
#  map_beta_div(Input_Image_File, Output_Dir, window_size, nb_partitions=nb_partitions,
#               nbCPU=nbCPU, MaxRAM=MaxRAM, nbclusters = nbclusters)

## ----alpha and beta diversity indices from vector layer------------------
#  # location of the spectral species raster needed for validation
#  TypePCA     = 'SPCA'
#  Dir.Raster  = file.path(Output.Dir,basename(Input.Image.File),TypePCA,'SpectralSpecies')
#  Name.Raster = 'SpectralSpecies'
#  Path.Raster = file.path(Dir.Raster,Name.Raster)
#  
#  # location of the directory where shapefiles used for validation are saved
#  vect        = system.file('extdata', 'VECTOR', package = 'biodivMapR')
#  Shannon.All = list() # ??
#  
#  # list vector data
#  Path.Vector         = list_shp(vect)
#  Name.Vector         = tools::file_path_sans_ext(basename(Path.Vector))
#  
#  # get alpha and beta diversity indicators corresponding to shapefiles
#  Biodiv.Indicators           = diversity_from_plots(Raster = Path.Raster, Plots = Path.Vector,NbClusters = nbclusters)
#  # if no name
#  Biodiv.Indicators$Name.Plot = seq(1,length(Biodiv.Indicators$Shannon[[1]]),by = 1)
#  Shannon.RS                  = c(Biodiv.Indicators$Shannon)[[1]]

## ----Write validation----------------------------------------------------
#  # write RS indicators
#  ####################################################
#  # write indicators for alpha diversity
#  Path.Results = file.path(Output.Dir, basename(Input.Image.File), TypePCA, 'VALIDATION')
#  dir.create(Path.Results, showWarnings = FALSE, recursive = TRUE)
#  ShannonIndexFile <- file.path(Path.Results, "ShannonIndex.tab")
#  write.table(Shannon.RS, file = ShannonIndexFile, sep = "\t", dec = ".", na = " ",
#              row.names = Biodiv.Indicators$Name.Plot, col.names= F, quote=FALSE)
#  
#  Results =  data.frame(Name.Vector, Biodiv.Indicators$Richness, Biodiv.Indicators$Fisher,                                Biodiv.Indicators$Shannon, Biodiv.Indicators$Simpson)
#  names(Results)  = c("ID_Plot", "Species_Richness", "Fisher", "Shannon", "Simpson")
#  write.table(Results, file = paste(Path.Results,"AlphaDiversity.tab",sep=''), sep="\t", dec=".",               na=" ", row.names = F, col.names= T,quote=FALSE)
#  
#  # write indicators for beta diversity
#  BC_mean = Biodiv.Indicators$BCdiss
#  colnames(BC_mean) = rownames(BC_mean) = Biodiv.Indicators$Name.Plot
#  write.table(BC_mean, file = paste(Path.Results,"BrayCurtis.csv",sep=''), sep="\t", dec=".", na=" ", row.names = F, col.names= T,quote=FALSE)
#  

## ----PCoA on Field Plots-------------------------------------------------
#  # apply ordination using PCoA (same as done for map_beta_div)
#  library(labdsv)
#  MatBCdist = as.dist(BC_mean, diag = FALSE, upper = FALSE)
#  BetaPCO   = pco(MatBCdist, k = 3)
#  

## ----plot PCoA & Shannon-------------------------------------------------
#  # very uglily assign vegetation type to polygons in shapefiles
#  nbSamples = c(6,4,7,7)
#  vg        = c('Forest high diversity', 'Forest low diversity', 'Forest medium diversity', 'low vegetation')
#  Type_Vegetation = c()
#  for (i in 1: length(nbSamples)){
#    for (j in 1:nbSamples[i]){
#      Type_Vegetation = c(Type_Vegetation,vg[i])
#    }
#  }
#  
#  # create data frame including alpha and beta diversity
#  library(ggplot2)
#  Results     =  data.frame('vgtype'=Type_Vegetation,'pco1'= BetaPCO$points[,1],'pco2'= BetaPCO$points[,2],'pco3' = BetaPCO$points[,3],'shannon'=Shannon.RS)
#  
#  # plot field data in the PCoA space, with size corresponding to shannon index
#  ggplot(Results, aes(x=pco1, y=pco2, color=vgtype,size=shannon)) +
#    geom_point(alpha=0.6) +
#    scale_color_manual(values=c("#e6140a", "#e6d214", "#e68214", "#145ae6"))
#  filename = file.path(Path.Results,'BetaDiversity_PcoA1_vs_PcoA2.png')
#  ggsave(filename, plot = last_plot(), device = 'png', path = NULL,
#         scale = 1, width = 6, height = 4, units = "in",
#         dpi = 600, limitsize = TRUE)
#  
#  
#  ggplot(Results, aes(x=pco1, y=pco3, color=vgtype,size=shannon)) +
#    geom_point(alpha=0.6) +
#    scale_color_manual(values=c("#e6140a", "#e6d214", "#e68214", "#145ae6"))
#  filename = file.path(Path.Results,'BetaDiversity_PcoA1_vs_PcoA3.png')
#  ggsave(filename, plot = last_plot(), device = 'png', path = NULL,
#         scale = 1, width = 6, height = 4, units = "in",
#         dpi = 600, limitsize = TRUE)
#  
#  ggplot(Results, aes(x=pco2, y=pco3, color=vgtype,size=shannon)) +
#    geom_point(alpha=0.6) +
#    scale_color_manual(values=c("#e6140a", "#e6d214", "#e68214", "#145ae6"))
#  filename = file.path(Path.Results,'BetaDiversity_PcoA2_vs_PcoA3.png')
#  ggsave(filename, plot = last_plot(), device = 'png', path = NULL,
#         scale = 1, width = 6, height = 4, units = "in",
#         dpi = 600, limitsize = TRUE)
#