Search code examples
rfor-looprastershapefile

How can I create a for loop function with multiple shapefile datasets?


I have created this script to compute one output from two constants rc_Mod2000_LC, y=S01_E031_FC.

# Extract LUC in Mod2000
LU_S01_E031_FC_Mod2000 <- extract(x=rc_Mod2000_LC, y=S01_E031_FC)
maj <- function(x){
  y <- as.numeric(names(which.max(table(x))))
  return(y)
}
LU_S01_E031_FC_Mod2000 <- extract(x=rc_Mod2000_LC, y=S01_E031_FC, fun=maj)
colnames (LU_S01_E031_FC_Mod2000) <- c("LU_2000")

Now, I want to create a for loop function where rc_Mod2000_LC is a constant (raster) but where S01_E031_FC is a variable (shapefile). rc_Mod2000_LC is a global land cover map and S01_E031_FC is a tile. I have this list of 61 shapefiles.

[[1]]
class       : SpatialPolygonsDataFrame 
features    : 16 
extent      : 30.95493, 31.02964, -1.040257, -0.9624111  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
variables   : 16
names       :   ID, Area_ORG, LU_1990, LU_2000, CHLU_90_00, LU_2005, CHLU_00_05,     Tile,       UNIQ_ID,      AREA, D_90_00, D_00_05, Sour_90_00, Sour_00_05, Conf_90_00, ... 
min values  : 1139,    10.44,      11,      11,       1111,       0,       1112, S01_E031, S01_E031_1139,   90260.9,     100,     200,    GLSE_00,    GLSE_05,          3, ... 
max values  :  935,     8.97,      11,      12,       1112,       0,       1230, S01_E031,  S01_E031_935, 1029736.2,     520,     520,    GLSE_00,    GLSE_05,          3, ... 

[[2]]
class       : SpatialPolygonsDataFrame 
features    : 2 
extent      : 30.95484, 30.99854, -2.022452, -1.971676  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
variables   : 16
names       :  ID, Area_ORG, LU_1990, LU_2000, CHLU_90_00, LU_2005, CHLU_00_05,     Tile,      UNIQ_ID,       AREA, D_90_00, D_00_05, Sour_90_00, Sour_00_05, Conf_90_00, ... 
min values  : 466,     0.06,      11,      11,       1111,       0,       1130, S02_E031, S02_E031_466,   603.7435,     100,     520,    GLSE_00,      GE_06,          3, ... 
max values  : 975,     3.15,      11,      11,       1111,       0,       1130, S02_E031, S02_E031_975, 31698.1260,     100,     520,    GLSE_00,    GLSE_05,          3, ... 

Therefore I want x=rc_Mod2000_LC and y=list of shapefiles. Hope somebody can help me out and let me know if you need some clarifications.


Solution

  • I would do something like this , First you extract function:

    extract_shape(y=shape,x=rc_Mod2000_LC){
      S01_E031_FC<- readShapePoly(shape)
      proj4string(S01_E031_FC) <- "+proj=longlat +datum=WGS84"
      extract(y=S01_E031_FC, fun=maj)
    }
    

    Then you loop through you shapes names:

    lapply(list_shapes,extract_shape)
    

    Where :

    list_shapes <- c("S01_E031_FC","S01_E032_FC","S01_E033_FC",...)
    

    Or better if your shapes names have a certain pattern you can use :

    list_shapes <-  ls(pattern='S1_.*_FC')