Search code examples
rgpsgisspatialr-sf

Displace GPS coordinates of sf point features by vector of offsets in R


Okay, so I'm building on a previous question. I've got an original set of data, and I'm following a specific protocol (that I need to follow for IRB reasons) to offset these randomly.

Each point has it's own randomly generated offset. For reference, urban points are offset by 0-2km, rural are offset by 0-5km with randomly selected 1% offset by 0-10km.

Here is my code so far, I'm using a projected coordinate reference system

displaceData <- function(sfDataset) {
  
  dataset.sf <- sfDataset
  
  # 2. Generate a random direction by generating angle between 0 and 360, and converting the angle from degrees to radians.
  dataset.sf$angle_rad <- deg2rad(runif(nrow(dataset.sf), 0, 360))
  
  # 3. Generate a random distance in meters of 0-5,000 meters for Rural points with 1% of rural points being given 0-10,000 meter distance.
  # get number of rural clusters
  nRuralClusters <- sum(dataset.sf$URBAN_RURA == "R")
  #   Split urban and rural
  dataset_R.sf <- subset(dataset.sf, URBAN_RURA == "R")
  dataset_U.sf <- subset(dataset.sf, URBAN_RURA == "U")
  #   For rural, assign 5000 meter displacement with 1% randomly assigned up to 10000 meter displacement
  dataset_R.sf$m_displaced <- ifelse({runif(nRuralClusters) |> rank(ties.method = "random") <= floor(nRuralClusters*0.01)}, 10000, 5000)
  #   For urban, assign 2000 meter displacement
  dataset_U.sf$m_displaced <- 2000
  #   Combine them back together
  dataset_C.sf <- rbind(dataset_R.sf, dataset_U.sf)
  
  dataset_C.sf$random_meters <- runif(nrow(dataset_C.sf), min = 0, max = dataset_C.sf$m_displaced)
  
  
  # 4. Generate the offset by applying trigonometry formulas (law of cosines) using the distance as the hypotenuse and the radians calculated in step 2.
  dataset_C.sf$xOffset <- sin(dataset_C.sf$angle_rad)*dataset_C.sf$random_meters
  dataset_C.sf$yOffset <- cos(dataset_C.sf$angle_rad)*dataset_C.sf$random_meters
  
  # 5. Add the offset to the original coordinate (in meters) to return the displaced coordinates.
    st_geometry(dataset_C.sf) <- st_geometry(dataset_C.sf) + cbind(dataset_C.sf$xOffset, dataset_C.sf$yOffset)
  
  #st_geometry(dataset_C.sf) <- st_sfc(st_point(c(st_coordinates(dataset_C.sf$geometry)[,1] + dataset_C.sf$xOffset, st_coordinates(dataset_C.sf$geometry)[,2] + dataset_C.sf$yOffset))) 
  
  return(dataset_C.sf)
}

Step 5 is where I'm getting hung up. I'm having trouble shifting the coordinates by the displacements.


Solution

  • I answered my own question after a while.

    
    displaceData <- function(sfDataset) {
      
      dataset.sf <- htDHS_t2.sf
      # dataset.sf <- sfDataset
      
      # 2. Generate a random direction by generating angle between 0 and 360, and converting the angle from degrees to radians.
      dataset.sf$angle_rad <- deg2rad(runif(nrow(dataset.sf), 0, 360))
      
      # 3. Generate a random distance in meters of 0-5,000 meters for Rural points with 1% of rural points being given 0-10,000 meter distance.
      # get number of rural clusters
      nRuralClusters <- sum(dataset.sf$URBAN_RURA == "R")
      #   Split urban and rural
      dataset_R.sf <- subset(dataset.sf, URBAN_RURA == "R")
      dataset_U.sf <- subset(dataset.sf, URBAN_RURA == "U")
      #   For rural, assign 5000 meter displacement with 1% randomly assigned up to 10000 meter displacement
      dataset_R.sf$m_displaced <- ifelse({runif(nRuralClusters) |> rank(ties.method = "random") <= floor(nRuralClusters*0.01)}, 10000, 5000)
      #   For urban, assign 2000 meter displacement
      dataset_U.sf$m_displaced <- 2000
      #   Combine them back together
      dataset_C.sf <- rbind(dataset_R.sf, dataset_U.sf)
      
      dataset_C.sf$random_meters <- runif(nrow(dataset_C.sf), min = 0, max = dataset_C.sf$m_displaced)
      
      
      # 4. Generate the offset by applying trigonometry formulas (law of cosines) using the distance as the hypotenuse and the radians calculated in step 2.
      dataset_C.sf$xOffset <- sin(dataset_C.sf$angle_rad)*dataset_C.sf$random_meters
      dataset_C.sf$yOffset <- cos(dataset_C.sf$angle_rad)*dataset_C.sf$random_meters
      
      # 5. Add the offset to the original coordinate (in meters) to return the displaced coordinates.
      dataset_C.sf$newX <- st_coordinates(dataset_C.sf)[,1] + dataset_C.sf$xOffset
      dataset_C.sf$newY <- st_coordinates(dataset_C.sf)[,2] + dataset_C.sf$yOffset
      
      # Remove geometry
      dataset_C <- st_set_geometry(dataset_C.sf, NULL)
      
      dataset_f.sf <- st_as_sf(dataset_C, coords = c("newX", "newY"), crs = st_crs(32618))
      
      #st_geometry(dataset_C.sf) <- st_sfc(st_point(c(st_coordinates(dataset_C.sf$geometry)[,1] + dataset_C.sf$xOffset, st_coordinates(dataset_C.sf$geometry)[,2] + dataset_C.sf$yOffset))) 
      
      return(dataset_f.sf)
    }
    

    Basically, I removed the geometry, and created a new geometry object using the updated points. Make sure to use the original CRS!