I am new to working with spatial data in R and I have run into a few issues, that I was hoping the community could help me with. Here is the situation: I generated a raster file of a country with a range of probabilities. Within that country, I have state-level data - # of farms. My goal is to randomly seed a certain number of points (farms) (# of farms is stored in CSV file, “Table2”) in each state of the country, based on the probability, indicated in the raster file. As such, I clipped my raster probability file up into 16 different raster files (the number of states in this country).
I have my code written, for what I describe above. But this code needs to be performed over each of the 16 layers - the ‘ncells’ and “TotalFarm’ value needs to be updated for each file.
install.packages("raster") install.packages("sp") install.packages("rgdal") setwd("~XXXXXX”) r <- raster("Suitability_Prob_v30.tif") #Reading in 16 rasters r_31 <- raster("Suitability_Prob_v31.tif") r_32 <- raster("Suitability_Prob_v32.tif") r_33 <- raster("Suitability_Prob_v33.tif") r_34 <- raster("Suitability_Prob_v34.tif") r_35 <- raster("Suitability_Prob_v35.tif") r_36 <- raster("Suitability_Prob_v36.tif") r_37 <- raster("Suitability_Prob_v37.tif") r_38 <- raster("Suitability_Prob_v38.tif") r_39 <- raster("Suitability_Prob_v39.tif") r_310 <- raster("Suitability_Prob_v310.tif") r_311 <- raster("Suitability_Prob_v311.tif") r_312 <- raster("Suitability_Prob_v31.tif") r_313 <- raster("Suitability_Prob_v313.tif") r_314 <- raster("Suitability_Prob_v314.tif") r_315 <- raster("Suitability_Prob_v315.tif") TABLE2 <- read.csv(file=“XXXXX”, header=TRUE, sep=",") #Read in cvs that has farm #'s
Here performing ONLY on the first of the 16 raster files Set NAs to 0
r[==-3.402823e+38] <- NA r[is.na(r)]=0
Calculating number of cells
Calling in the # of farms for this state, stored in CSV
For point selection, set number of cells (ncells) and the number of points to be selected, (specific to each state) Weighted by the value in the cells:
ptscell = sample(ncells, TotalFarm, prob=r, replace=TRUE)
Distribute the points throughout the grid cells Get the cell half-width:
hs = res(r)/2
Now find the center of those points:
centers = xyFromCell(r,ptscell)
Generate random uniform points in the cell by using the center and the half-width/height
pts = cbind(runif(nrow(centers),centers[,1]-hs,centers[,1]+hs), runif(nrow(centers),centers[,2]-hs,centers[,2]+hs))
Print and finish:
plot(r) points(pts) print(pts) write.csv(pts, file = "Points_5.txt",row.names=FALSE)
My attempts to apply this code over the series of 16 states has failed, in part, because I cannot get my 16 raster files into a “raster object” that allows me to use the rasters individually with those changes. Specifically, I made a mosaic image out of the 16 rasters , however the changes that I made on the mosaic were not made on the original raster files (only on the new mosaic object), so when I need to go back and seed the files individually with the number of farms specific to each state, the modifications that I do to the mosaic, i.e. NA’s, are not there.
All of the tutorials that I have seen perform loops on rasters by building a brick or stack. But since these are not bands which have similar extents, rather, they are neighboring states, this approach is not appropriate here.
So, how do I perform a loop on rasters, that are in the “raster” class or otherwise apply the code written above, to each raster automatically?