[GIS] How to calculate mean for multiple raster files and multiple variables using R

averagerrasterstack

Folder with multiple raster filesI have a folder with multiple raster files (.asc) which contains the biomass of different fish species (e.g., Anglerfish, Bathydemersal etc) for different years (e.g., 00012, 00024, 00036 etc)
I want to calculate the mean per species using R
So far i have done the following:

files_path <- "~/new/asc/" #path where my files are    
all_files <- list.files(files_path,
                            full.names = TRUE,
                            pattern = ".asc$") #take all the ascii files in    
files_stack <- stack(all_files) #stack them together    

Here comes the problem: mean_ <- calc(files_stack, mean) code makes one mean for the entire files_stack. Mean for all the species in my stack

I want to have mean per each species. Any tips on how to do it?

Best Answer

Once you have a stack of rasters loaded, you can operate on groups of them using stackApply. You have to create a vector of the groups.

In your case the grouping is the species. Extract the species from the file name using some string processing. For example if your file paths are like this:

> all_files
 [1] "./EcospaceMapBiomass-Anglerfish-00012.asc"   
[etc]
 [7] "./EcospaceMapBiomass-Anglerfish-00084.asc"   
 [8] "./EcospaceMapBiomass-Bethydemersal-00012.asc"
 [9] "./EcospaceMapBiomass-Bethydemersal-00024.asc"
[10] "./EcospaceMapBiomass-Bethydemersal-00036.asc"
[etc]      
[18] "./EcospaceMapBiomass-Shark-00048.asc"        
[19] "./EcospaceMapBiomass-Shark-00060.asc"        
[20] "./EcospaceMapBiomass-Shark-00072.asc"        
[21] "./EcospaceMapBiomass-Shark-00084.asc"  

Then you can get the species by removing everything up to "mass-" with nothing, and then everything from the remaining "-" to the end with nothing. This will probably fail if you have a "-" in a species name. Check that. So this:

species = gsub("-.*.asc","",gsub(".*mass-","",all_files))

nets you:

> species
 [1] "Anglerfish"    "Anglerfish"    "Anglerfish"    "Anglerfish"   
 [5] "Anglerfish"    "Anglerfish"    "Anglerfish"    "Bethydemersal"
 [9] "Bethydemersal" "Bethydemersal" "Bethydemersal" "Bethydemersal"
[13] "Bethydemersal" "Bethydemersal" "Shark"         "Shark"        
[17] "Shark"         "Shark"         "Shark"         "Shark"        
[21] "Shark"    

You could also do this with strsplit or one of a number of other ways. However you do it, end up with a vector of species that match each item in your list of files.

Then stackApply.

> files_stack <- stack(all_files) 
> means = stackApply(files_stack, indices=species, fun=mean)
> means
class      : RasterBrick 
dimensions : 4, 4, 16, 3  (nrow, ncol, ncell, nlayers)
resolution : 0.25, 0.25  (x, y)
extent     : 0, 1, 0, 1  (xmin, xmax, ymin, ymax)
crs        : NA 
source     : memory
names      : index_Anglerfish, index_Bethydemersal, index_Shark 
min values :      0.006978279,         0.006978279, 0.006978279 
max values :        0.9034359,           0.9034359,   0.9034359 

means is now a brick with 3 layers, each one a map of the mean for a species, and the layer names give you the species. These are all identical because all my inputs were identical. Test with better simulated data to make sure this is actually correct, or make sure the answer with real data looks reasonable.