Ответ 1

Я придумал подход, который был бы применим по сравнению с многочисленными вызовами API.

Идея состоит в том, чтобы найти места, куда вы можете добраться за определенное время (посмотрите эту ветку) Трафик можно смоделировать, изменив время с утра до вечера. Вы получите перекрывающуюся область, в которую можно попасть из обоих мест.

Затем вы можете использовать ответ Николаса и нанести на карту некоторые точки в пределах этой перекрывающейся области и нарисовать карту тепла для пунктов назначения, которые у вас есть. Таким образом, у вас будет меньше площади (точек) для покрытия, и, следовательно, вы будете делать гораздо меньше вызовов API (не забудьте использовать соответствующее время для этого вопроса).

Ниже я попытался продемонстрировать, что я имею в виду под этим, и привести вас к выводу, что вы можете сделать сетку, упомянутую в другом ответе, чтобы сделать вашу оценку более надежной.

Это показывает, как отобразить пересеченную область.

library(httr)
library(googleway)
library(jsonlite)

appId <- "Travel.Time.ID"
apiKey <- "Travel.Time.API"
mapKey <- "Google.Map.ID"




locationK <- c(40, -73) #K
locationM <- c(40, -74) #M

CommuteTimeK <- (3 / 4) * 60 * 60
CommuteTimeM <- (0.55) * 60 * 60

url <- "http://api.traveltimeapp.com/v4/time-map"

requestBodyK <- paste0('{ 
                      "departure_searches" : [ 
                      {"id" : "test", 
                      "coords": {"lat":', locationK[1], ', "lng":', locationK[2],' }, 
                      "transportation" : {"type" : "public_transport"} ,
                      "travel_time" : ', CommuteTimeK, ',
                      "departure_time" : "2018-06-27T13:00:00z"
                      } 
                      ] 
                      }')


requestBodyM <- paste0('{ 
                      "departure_searches" : [ 
                      {"id" : "test", 
                      "coords": {"lat":', locationM[1], ', "lng":', locationM[2],' }, 
                      "transportation" : {"type" : "driving"} ,
                      "travel_time" : ', CommuteTimeM, ',
                      "departure_time" : "2018-06-27T13:00:00z"
                      } 
                      ] 
                      }')


resKi <- httr::POST(url = url,
                  httr::add_headers('Content-Type' = 'application/json'),
                  httr::add_headers('Accept' = 'application/json'),
                  httr::add_headers('X-Application-Id' = appId),
                  httr::add_headers('X-Api-Key' = apiKey),
                  body = requestBodyK,
                  encode = "json")


resMi <- httr::POST(url = url,
                   httr::add_headers('Content-Type' = 'application/json'),
                   httr::add_headers('Accept' = 'application/json'),
                   httr::add_headers('X-Application-Id' = appId),
                   httr::add_headers('X-Api-Key' = apiKey),
                   body = requestBodyM,
                   encode = "json")

resK <- jsonlite::fromJSON(as.character(resKi))
resM <- jsonlite::fromJSON(as.character(resMi))

plK <- lapply(resK$results$shapes[[1]]$shell, function(x){
  googleway::encode_pl(lat = x[['lat']], lon = x[['lng']])
})

plM <- lapply(resM$results$shapes[[1]]$shell, function(x){
  googleway::encode_pl(lat = x[['lat']], lon = x[['lng']])
})

dfK <- data.frame(polyline = unlist(plK))
dfM <- data.frame(polyline = unlist(plM))

df_markerK <- data.frame(lat = locationK[1], lon = locationK[2], colour = "#green")
df_markerM <- data.frame(lat = locationM[1], lon = locationM[2], colour = "#lavender")

iconK <- "red"
df_markerK$icon <- iconK

iconM <- "blue"
df_markerM$icon <- iconM


google_map(key = mapKey) %>%
  add_markers(data = df_markerK,
              lat = "lat", lon = "lon",colour = "icon",
              mouse_over = "K_K") %>%
  add_markers(data = df_markerM, 
              lat = "lat", lon = "lon", colour = "icon",
              mouse_over = "M_M") %>%
  add_polygons(data = dfM, polyline = "polyline", stroke_colour = '#461B7E',
               fill_colour = '#461B7E', fill_opacity = 0.6) %>% 
  add_polygons(data = dfK, polyline = "polyline", 
               stroke_colour = '#F70D1A',
               fill_colour = '#FF2400', fill_opacity = 0.4)

Вы можете извлечь пересеченную область следующим образом:

install.packages(c("rgdal", "sp", "raster","rgeos","maptools"))
library(rgdal)
library(sp)
library(raster)
library(rgeos)
library(maptools)

Kdata <- resK$results$shapes[[1]]$shell
Mdata <- resM$results$shapes[[1]]$shell


xyfunc <- function(mydf) {
  xy <- mydf[,c(2,1)]
  return(xy)
}

spdf <- function(xy, mydf) {sp::SpatialPointsDataFrame(coords = xy, data = mydf,
                                                       proj4string = CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0"))}


for (i in (1:length(Kdata))) {Kdata[[i]] <- xyfunc(Kdata[[i]])}

for (i in (1:length(Mdata))) {Mdata[[i]] <- xyfunc(Mdata[[i]])}


Kshp <- list()
for (i in (1:length(Kdata))) {Kshp[i] <- spdf(Kdata[[i]],Kdata[[i]])}

Mshp <- list()
for (i in (1:length(Mdata))) {Mshp[i] <- spdf(Mdata[[i]],Mdata[[i]])}

Kbind <- do.call(bind, Kshp) 
Mbind <- do.call(bind, Mshp) 
#plot(Kbind)
#plot(Mbind)


x <- intersect(Kbind,Mbind)
#plot(x)

xdf <- data.frame(x)
head(xdf)
#         lng      lat     lng.1    lat.1 optional
# 1 -74.23374 40.77234 -74.23374 40.77234     TRUE
# 2 -74.23329 40.77279 -74.23329 40.77279     TRUE
# 3 -74.23150 40.77279 -74.23150 40.77279     TRUE
# 4 -74.23105 40.77234 -74.23105 40.77234     TRUE
# 5 -74.23239 40.77099 -74.23239 40.77099     TRUE
# 6 -74.23419 40.77099 -74.23419 40.77099     TRUE


xdf$icon <- "https://i.stack.imgur.com/z7NnE.png"
google_map(key = mapKey, location = c(mean(latmax,latmin), mean(lngmax,lngmin)), zoom = 8) %>% 
     add_markers(data = xdf, lat = "lat", lon = "lng", marker_icon = "icon")

Это просто иллюстрация пересекаемой области.

Теперь вы можете получить координаты из xdf и построить свою сетку вокруг этих точек, чтобы в итоге получить тепловую карту. Чтобы уважать другого пользователя, который придумал эту идею/ответ, я не включаю его в свой, а просто ссылаюсь на него.

Николас Веласкес - Получение матрицы происхождения-назначения между сеткой (грубо) одинаково удаленных точек

Ответ 2

Этот ответ основан на получении матрицы отправления-назначения между сеткой (примерно) одинаково удаленных точек. Это ресурсоемкая операция не только потому, что она требует большого числа API-вызовов для картографических сервисов, но и потому, что серверы должны вычислять матрицу для каждого вызова. Количество необходимых вызовов растет экспоненциально по количеству точек в сетке.

Чтобы решить эту проблему, я бы посоветовал вам запустить на локальном компьютере или на локальном сервере сервер сопоставления. Project OSRM предлагает относительно простое бесплатное решение с открытым исходным кодом, позволяющее запускать сервер OpenStreetMap в докер Linux (https://github.com/Project-OSRM/osrm-backend). Наличие собственного локального сервера сопоставления позволит вам совершать столько вызовов API, сколько пожелаете. Пакет R osrm позволяет вам взаимодействовать с apis OpenStreetMaps. в том числе размещены на локальном сервере.

library(raster) # Optional
library(sp)
library(ggmap)
library(tidyverse)
library(osrm)
devtools::install_github("cmartin/ggConvexHull") # Needed to quickly draw the contours
library(ggConvexHull)

Я создаю сетку из 96 примерно одинаково удаленных точек вокруг Брюсселя (Бельгия). Эта сетка не учитывает кривизну земли, которая незначительна на уровне городских расстояний.

Для удобства я использую растровый пакет для загрузки ShapeFile из Бельгии и извлечения узлов для города Брюссель.

  BE <- raster::getData("GADM", country = "BEL", level = 1)
  Bruxelles <- BE[BE$NAME_1 == "Bruxelles", ]

  df_grid <- makegrid(Bruxelles, cellsize = 0.02) %>% 
        SpatialPoints() %>% 
        as.data.frame() %>% ## I convert the SpatialPoints object into a simple data.frame
        rownames_to_column() %>% ## create a unique id for each point in the data.frame
        rename(id = rowname, lat = x2, lon = x1) # rename variables of the data.frame with more explanatory names.

 options(osrm.server = "http://127.0.0.1:5000/") ## I point osrm.server to the OpenStreet docker running in my Linux machine. Do not run this if you are getting your data from OpenStreet public servers.

 Distance_Tables <- osrmTable(loc = df_grid)  ## I obtain a list with distances (Origin Destination Matrix in minutes, origins and destinations)

 OD_Matrix <- Distance_Tables$durations %>% ## Subset the previous list and 
   as_data_frame() %>%  ## ...convert the Origin Destination Matrix into a tibble
   rownames_to_column() %>% 
   rename(origin_id = rowname) %>% ## make sure we have an id column for the OD tibble
   gather(key = destination_id, value = distance_time, -origin_id) %>% # transform the tibble into long/tidy format
   left_join(df_grid, by = c("origin_id" = "id")) %>% 
   rename(origin_lon = lon, origin_lat = lat) %>% ## set origin coordinates
   left_join(df_grid, by = c("destination_id" = "id")) %>% 
   rename(destination_lat = lat, destination_lon = lon) ## set destination coordinates

 ## Obtain a nice looking road map of Brussels

 Brux_map <- get_map(location = "bruxelles, belgique", 
                     zoom = 11, 
                     source = "google", 
                     maptype = "roadmap")

 ggmap(Brux_map) + 
   geom_point(aes(x = origin_lon, y = origin_lat), 
         data = OD_Matrix %>% 
                filter(destination_id == 42), ## Here I selected point_id 42 as the desired target, just because it is not far from the City Center.
                size = 0.5) + 
   geom_point(aes(x = origin_lon, y = origin_lat), 
        data = OD_Matrix %>% 
        filter(destination_id == 42, origin_id == 42),
          shape = 5, size = 3) +  ## Draw a diamond around point_id 42                                      
   geom_convexhull(alpha = 0.2, 
         fill = "blue", 
         colour = "blue",
         data = OD_Matrix %>% 
                filter(destination_id == 42, 
                       distance_time <= 8), ## Countour marking a distance of up to 8 minutes
         aes(x = origin_lon, y = origin_lat)) + 
   geom_convexhull(alpha = 0.2, 
         fill = "red",
         colour = "red",
         data = OD_Matrix %>% 
         filter(destination_id == 42, 
                distance_time <= 15), ## Countour marking a distance of up to 16 minutes
         aes(x = origin_lon, y = origin_lat))

Результаты

Синий контур обозначает расстояние до центра города до 8 минут. Красный контур отображает расстояния до 15 минут.

Я надеюсь, что это поможет вам получить ваши обратные изохроны.