Function for fast spatial thinning of geographic coordinates

thin_coordinates.R

thin_coordinates = function(coords, threshold){
  # Function for faster spatial thinning of coordinates
  # coords: geographic coordinates in lon/lat format
  # threshold: minimum distance between points in meters
  
  coords_dist = geosphere::distm(coords, fun = distHaversine)   # Great circle distance between points
  coords_dist[upper.tri(coords_dist, diag = T)] = NA            # Triangular distance matrix
  coords_flagged = which(coords_dist < threshold, arr.ind=TRUE) # Two-column matrix of coordinate-pairs that are below threshold
  pts_remove = c()
  while(nrow(coords_flagged) > 0){
    pts_count = sort(table(c(coords_flagged)))                  # Count number of occurrences of points in flagged coords
    pts_count_max = pts_count[pts_count == max(pts_count)]      # Points with most neighbors below threshold
    pt_target = as.integer(sample(names(pts_count_max), 1))     # Choose one randomly
    pt_ind = which(coords_flagged[,1] == pt_target | coords_flagged[,2] == pt_target) # Find all pairs that involve target point
    coords_flagged = coords_flagged[-pt_ind,,drop = F]          # Remove those pairs from coords_flagged
    pts_remove = c(pts_remove, pt_target) 
  }
  return(coords[-pts_remove,])
}