I want to generate 5 pairs of random points per each stratum. I have 4 strata that represent different rock types. For each rock type i want to generate 5 random points and to each random point I want to assign its respective pair at 10 meters distance (but it should still belong to the same stratum/category of rocks). So, altogether there should be 20 original random points and 20 additional (paired ones).

I have tried using sampling::strata (R package is called sampling and to get random points stratified per category the function strata is appropriate to use). See the code below. My sf object (data) is called selected_rows and the column with categories/strata is Rock_name. The size is 5 points per stratum.

units <- sampling::strata(
selected_rows, stratanames = "Rock_name", size = n_s[stratum], method = "srswr")

points = getdata(selected_rows, units) %>% st_as_sf()

My issue is now that I do not get coordinates of the random points rather the geometry is multi-polygon. I suppose from the strata. I would like to have coordinates of each of the 20 generated points so I can generate their respective pairs at 10 meters distance with random azimuth or st_buffer.

Does anyone know if it is possible to get the coordinates of each of the randomly generated points per stratum?

I was following the instructions from Dick Brus new book on spatial sampling. See the screenshot below. enter image description here


2 Answers 2


You can do that with terra like this:

Example data

v <- vect(system.file("ex/lux.shp", package="terra"))


s <- spatSample(v, 5, method="random", strata="ID_2")

plot(v); points(s, col="red")

One thing you could do with polygon data is first create a random sample that is larger than you need, using sf::st_sample, then use the stratified.random function in the spatialEco package. All you need is a stratifying variable and the function will create a stratified subsample that is an spatial sf object.

There is also a subsample.distance function available in the same package that subsamples based on a conditional min or max distance.

It is looking like you may be using a raster as the data that the sample is drawn from. Here is an example (if you can process your raster in memory) of a distance-based stratified sample using subsample.distance and lapply (one could also use for).

Add packages and create some example data, "r" is the raster that the stratification is drawn from.


r <- rast(ext(485660.1887, 898660.1887, 4873781.3587, 5113781.3587), 
          resolution = 1000, crs="ESRI:54030")
  r[] <- round(runif(ncell(r), 1, 2),0)
    r <- app(focal(r, gaussian.kernel(sigma=1, s=33)), fun=round)

# available strata counts

Coerce the raster to an sf POINT object and assign "STRAT" as the column name that will define the stratifying values.

rsf <- st_as_sf(as.points(r))
  names(rsf)[1] <- "STRAT"

Now, we can iterate through the strata using lapply and create a distance based random sample (see help for subsample.distance to understand the various argument).

samp <- lapply(unique(rsf$STRAT), function(x) {
  subsample.distance(rsf[rsf$STRAT == x,], size = 10, d = 5000)   
  samp <- do.call("rbind", samp)
  plot(samp["STRAT"], pch=20, cex=2, add=TRUE)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.