The raster package for R provides a range of GIS-like functions for analysing spatial grid data. Together with package sp, and several other spatial analysis packages, R provide a quite comprehensive set of tools for manipulating and analysing spatial data.

I needed to randomly select some locations for field sampling, with inclusion probabilities based on values contained in a raster. The code below did the job very easily.

library(raster) #an example raster from the raster package f <- system.file("external/test.grd", package="raster") r<-raster(f) plot(r) #make a raster defining the desired inclusion probabilities #for the all locations available for sampling probrast<-raster(r) #inclusion probability for cells with value >=400 #will be 10 times that for cells with value <400 probrast[r>=400]<-10 probrast[r<400]<-1 #normalise the probability raster by dividing #by the sum of all inclusion weights: probrast<-probrast/sum(getValues(probrast), na.rm=T) #confirm sum of probabilities is one sum(getValues(probrast), na.rm=T) #plot the raster of inclusion probabilities plot(probrast, col=c(gray(0.7), gray(0.3))) #a function to select N points on a raster, with #inclusion probabilities defined by the raster values. probsel<-function(probrast, N){ x<-getValues(probrast) #set NA cells in raster to zero x[is.na(x)]<-0 samp<-sample(nrow(probrast)*ncol(probrast), size=N, prob=x) samprast<-raster(probrast) samprast[samp]<-1 #set value of sampled squares to 1 #convert to SpatialPoints points<-rasterToPoints(samprast, fun=function(x){x>0}) s points<-SpatialPoints(points) return(points) } #select 300 sites using the inclusion probabilities #defined in probrast samppoints<-probsel(probrast, 300) plot(probrast, col=c(gray(0.7), gray(0.3)), axes=F) plot(samppoints, add=T, pch=16, cex=0.8, col="red")

Here’s the result. Note the higher density of sampled points (red) within the parts of the raster with higher inclusion probability (dark grey).

Advertisements