How do I find the area under different classification of map using raster data?

Hello,
Is it possible to estimate the area under different regions of the map obtained below preferably in % of total area?

Thank you

suppressPackageStartupMessages({
  library(tmap)
  library(tidyverse)
  library(rgdal)
  library(rgeos)
  library(spatstat)  
  library(maptools)  
  library(raster)    
  library(gstat) 
  library(sp)    
  library(gsheet)
})
P_month=gsheet2tbl('https://docs.google.com/spreadsheets/d/1rXB9rG4aRX8Hzyc_O97QcAII9AsiZDlgQ-1JsnMf6ww/edit?usp=sharing')
coordinates(P_month) <- ~lon+lat
proj4string(P_month) <- CRS("+init=epsg:28992")
States=readRDS(url("https://biogeo.ucdavis.edu/data/gadm3.6/Rsp/gadm36_USA_1_sp.rds","rb"))
W  <- States %>% 
  subset(NAME_1%in%c("Alabama","Arkansas", "Florida", "Georgia", "Kentucky", "Louisiana", "Mississippi","Tennessee"))
P_month@bbox <- W@bbox
th  <-  as(dirichlet(as.ppp(P_month)), "SpatialPolygons")
proj4string(th) <- proj4string(P_month)
grd              <- as.data.frame(spsample(P_month, "regular", n=1000000))
names(grd)       <- c("X", "Y")
coordinates(grd) <- c("X", "Y")
gridded(grd)     <- TRUE 
fullgrid(grd)    <- TRUE  
proj4string(P_month) <- proj4string(P_month) 
proj4string(grd) <- proj4string(P_month)
P_month$X <- coordinates(P_month)[,1]
P_month$Y <- coordinates(P_month)[,2]
f.1 <- as.formula(rain ~ X + Y) 
var.smpl <- variogram(f.1, P_month, cloud = FALSE)
dat.fit  <- fit.variogram(var.smpl, vgm(c("Exp", "Sph", "Gau", "Mat")))


dat.krg <- krige( f.1, P_month, grd, dat.fit)
r.m <- mask(raster(dat.krg), W)
plot(var.smpl, dat.fit)

tm_shape(r.m) +
  tm_raster(n=5, midpoint = NA,
            title="Rain",
          legend.is.portrait = FALSE) +
  tm_layout(legend.position=c('left', 'bottom'),
            frame = TRUE,
            legend.outside = FALSE,
            outer.margins=0)

You can simply count the cells within a given interval. Below two ways you could go, directly using the raster object or transforming to sf_points. Be aware that you need to use either the grid-resolution res() to arrive at actual area or divide by total cell-count for relative area, as you have equal distance grid.

Hope this gets you going.
JW

#directly using raster                                                                                   
#reclassification of spatraster using the bands above                                                    
recalc  <- matrix(c(0, 1000, 1, 1000, 2000, 2, 2000, 3000, 3, 3000, 4000, 4, 4000, 5000, 5, 5000, 10000, 
6), ncol = 3, byrow = TRUE)                                                                              
rm2  <- reclassify(r.m, recalc)                                                                          
freq(rm2, useNA = 'no')  #use frequency function                                                         
#      value  count                                                                                      
# [1,]     1 238752                                                                                      
# [2,]     2 126221                                                                                      
# [3,]     3  50368                                                                                      
# [4,]     4  31189                                                                                      
# [5,]     5   8763                                                                                      
# [6,]     6   1631                                                                                      

Alternatively you can convert all grid cells to sf_points and use dplyr verbs for further tranformation.

# convert to vector representation and use dplyr verbs to work on the data                   
sto  <- stars::st_as_stars(r.m)     #library stars to transform raster                       
# now make each raster cell a single sf_point                                                
rm_sf_points  <- st_as_sf(sto, as_points = TRUE, merge = FALSE, long = TRUE)  #              
                                                                                             
#add column with similar pixel-values the point atrributes. Here following your intevals     
rm_sf_points  <- rm_sf_points  %>%                                                           
        mutate(band = case_when(                                                             
                                var1.pred < 1000 ~ 1L,                                       
                                var1.pred  >= 1000 & var1.pred < 2000 ~ 2L,                  
                                var1.pred >= 2000 & var1.pred < 3000 ~ 3L,                   
                                var1.pred  >= 3000 & var1.pred < 4000 ~ 4L,                  
                                var1.pred  >= 4000 & var1.pred < 5000 ~ 5L,                  
                                var1.pred  >= 5000 ~ 6L))                                    
                                                                                             
#tm_shape(rm_sf_points) +tm_dots(col = "band")    #check to see                              
#and do the aggregation by regular grouping and aggregating                                  
freq_tab  <- rm_sf_points  %>% st_drop_geometry()  %>%                                       
        group_by(band)  %>%                                                                  
        tally()                                                                              
freq_tab                                                                                     
# # A tibble: 6 × 2                                                                          
#    band      n                                                                             
#   <int>  <int>                                                                             
# 1     1 238752                                                                             
# 2     2 126221                                                                             
# 3     3  50368                                                                             
# 4     4  31189                                                                             
# 5     5   8763                                                                             
# 6     6   1631                                                                             
                                                                                             

This topic was automatically closed 21 days after the last reply. New replies are no longer allowed.

If you have a query related to it or one of the replies, start a new topic and refer back with a link.