How can I efficiently aggregate with time and distance windows with data.table?

user_2402465 0
Viewed 99

I'm working with GPS data in a simulation study, which means that all my minor inefficiencies are coming back to haunt me.

One of the bigger issues is implementation of a stop detection algorithm. The first step is the generation of candidate stops, where I identify coordinates that are close together in time and space.

Here's some toy data:

library(data.table)
gcdist <- function(lat1, lon1, lat2, lon2){
  sqrt((12430*abs(lat1 - lat2)/180)^2 +
         ((24901*abs(lon1 - lon2)/360) * cos((lat1+lat2)/2))^2) 
}

time <- rev(Sys.time() - seq(60, 30*60, 60))
lat <- cumsum(c(55, rnorm(29, 0, sd  = .003)))
lon <- cumsum(c(5.2, rnorm(29, 0, sd  = .003)))
id <- 1:30

dt <-data.table(id, time, time_dup = time, lon, lat)
dt[, time_window := time + 60*3]
setkey(dt, time, time_window)

What I'm currently doing uses foverlap from data.table. Data.table-friendly solutions are kind of a must because the entire data set approaches 8GB.

This is not efficient enough to run my simulation study on a decent timeline, and it's also not very precise. Because I only look within the three minute window and then stop, I have to find ways to handle combining them.

# Current setup
setkey(dt, time, time_window)
temp <- foverlaps(dt, dt, by.y = c("time", "time_window"), by.x = c("time", "time_dup"), type = "within")
temp[, dist := gcdist(lat, lon, i.lat, i.lon)]
temp[dist < .5, within_stop  := TRUE]
temp[, candidate_stop := all(dist < .5), id]
setkey(temp, id)
setkey(dt, id)
dt[temp, candidate_stop := i.candidate_stop]
setkey(temp, i.id)
dt[temp, within_stop := i.within_stop]

This is closer to what I ought to be doing, but it's just too cumbersome for 8gb of data.

current_stop <- dt[1, ]  

for (i in seq_len(nrow(dt))) {
  
  dist <- gcdist(current_stop[, lat], current_stop[, lon],
       dt[i, lat], dt[i, lon])
  if(abs(dist) > .5) current_stop <- dt[i]
  dt[i:nrow(dt), stop_id := current_stop[, id]]
}

I have a hunch there's something for me in data.table's rolling joins or frollapply. I keep reading through the pages and toying with examples, but I can't quite make it work. I think I should be able to rollapply a distance windowing function until the distance is greater than a certain cutoff, then restart, but beats me if I can figure out how.

rolling-computation data.table gps r
1 Answers

I believe this should help, although there may be some issues with your within_stop variable.

dt[dt,
   on = .(time >= time, time_dup <= time_window),
   c("within_stop", "candidate_stop") := {
     within_stop = gcdist(lat, lon, i.lat, i.lon) < 0.5
     .(within_stop = first(within_stop), ## I am least sure about this. 
          candidate_stop = all(within_stop))
   },
   by = .EACHI]

dt

On my machine, this is about twice as fast on this dataset. I mention within_stop because it likely varies within each join group. That is, sometimes it is all TRUE while other times it could be a mixture. I could not determine what dt[dt, :=] defaults to if there is multiple matches although I would assume it is likely to be first.

To make it slightly quicker, you likely could use compiled code with rcpp in which you would implement gcdist() with short-circuiting so that it only returns a boolean. That way you would not have to allocate a vector of distances and a boolean vector of whether it is within your tolerance. This seems to work although there is very little performance gain for me with your current dataset.

bool gcdist_short_circuit(NumericVector lat, NumericVector lon,double i_lat, double i_lon) {

  bool out = TRUE;

  for (int i = 0; i < lat.size(); i++) {
    if (sqrt(pow(12430*abs(lat[i] - i_lat)/180, 2) + pow((24901*abs(lon[i] - i_lon)/360) * cos((lat[i]+i_lat)/2),2)) >= 0.5) {
       out = FALSE;
       break;
     }
   }
  return(out);
}
user_113550660
Related