Получение открытых данных Open Street. Map для оценки запечатанности парковых территорий

1.

Получение открытых
данных Open Street
Map для оценки
запечатанности
парковых территорий
Ярославцев А.М.

2.

Open Street Map
2
OpenStreetMap (дословно «открытая карта улиц»),
сокращённо OSM — некоммерческий вебкартографический проект по созданию силами
сообщества участников — пользователей Интернета
подробной свободной и бесплатной географической
карты мира. https://www.openstreetmap.org/.
Готовые для использования карты доступны онлайн
на самом сайте OpenStreetMap и через множество
других сайтов. С помощью функции «Экспорт»
можно получить HTML-код для вставки на любой
сайт с возможностью нанесения маркера или
непосредственные данные в формате OSM(вариант
2022-03-15

3.

Пакеты для работы с OSM
3
#load packages
library(tidyverse)
library(sf)
sf::sf_use_s2(FALSE)
library(osmdata)
library(ggmap)
library(leaflet)
library(RColorBrewer)
Нам понадобятся пакеты tidyverse и sf для общих и
пространственных манипуляций над данными
соответственно. sf::sf_use_s2(FALSE) - включает
настройку позволяющую проводить операции с
2022-03-15

4.

Загрузим полигон с
границами парка
park_sf =
read_sf('shapes/dubky.geojson')
plot(park_sf)
4
2022-03-15

5.

Предварительный
этап
get_overpass_url()
Исходно, пакет osmdata
## [1]
"https://overpassapi.de/api/interpreter
"
overpass_url =
"https://maps.mail.ru/
osm/tools/overpass/api
/interpreter"
set_overpass_url(overp
ass_url)
bbox =
st_bbox(park_sf$geomet
ry) %>%
matrix(ncol=2,nrow=2)
5
2022-03-15
получает данные не от
OSM, а от одного из
зеркал. Вы можете
установить в качестве
источника данных, другое
зеркало, с меньшей
нагрузкой и задержками.
Чтобы получить данные
интересующего вас
объекта, пакету osmdata
необходимо сообщить
границы области объекта
в виде матрицы

6.

Данные
OSM
от
GGmap
mad_map <- get_map(bbox,
mad_map <- get_map(bbox,
zoom = 17, source =
zoom = 16, source =
разнятся
по
доступности на
"stamen")
"stamen")
plot(mad_map)
plot(mad_map)
разных масштабах
6
2022-03-15

7.

Получение данных OSM
bbox
=
Для OSMDATA
данных через
через
пакет
st_bbox(park_sf$geomet
библиотеку OSMDATA мы
ry) %>%
matrix(ncol=2,nrow=2)
colnames(bbox) =
c("min","max")
rownames(bbox) =
c("x","y")
dubky_all = bbox %>%
opq(timeout = 900) %>%
add_osm_feature(key =
"natural",
value =
available_tags("natura
7
2022-03-15
l")) %>%
должны получить
координаты объекта в
виде матрицы, с
заданными именами строк
и столбцов. Далее эта
матрица используется в
функции opq, которая
обращается к серверу и
запрашивает данные из
категории “natural” -т.е о
естественных объектах,
чтобы получить все
данные из этой категории
мы используем команду
available_tags(“natural”),
после чего конвертируем
все в пространственный

8.

Получение данных OSM
ggplot()+
Полученные
данные легко
через
пакет
OSMDATA
geom_sf(data =
отобразить с помощью
dubky_all$osm_polygons
, aes())+
geom_sf(data =
dubky_all$osm_multipol
ygons,
aes(fill=natural))+
theme_bw()
8
2022-03-15
пакета ggplot и в
частности команды
geom_sf(), позволяющей
отображать в виде карт
содержимое объектов
класса sf. К сожалению
при автоматизированных
запросах, даже на
альтернативных серверах,
часто в качестве ответа
можно получить только
часть данных - как в
данном случае. Мы
получили данные только
по зеленым насаждениям.
Обойти эту ситуацию
можно только

9.

Получение данных OSM
напрямую
Сайт openstreetmaps.org позволяет получать полную
выгрузку, со всеми данными объекта, в разделе
export. Скачиваемый файл формата .osm
представляет собой данные структурированные в
формаие XML
9
2022-03-15

10.

Подготовка данных OSM к
map_polygons
=
К счастью получать
анализу
st_read("map.osm",
данные из файлов
layer =
'multipolygons', quiet
= TRUE)
map_lines =
st_read("map.osm",
layer = 'lines', quiet
= TRUE)
map_polygons =
map_polygons %>%
filter(!is.na(natural)
|
!
is.na(building))
map_lines = map_lines
%>%
filter(highway
==
10
2022-03-15
"footway")
формата osm могут
стандартные функции
чтения пространственных
данных из пакета sf st_read(). Главный нюанс
их использования
заключается в том, что
нельзя получить все типы
пространственной
геометрии из одного
файла и их можно
“вытянуть” только
последовательно в
различные объекты. При
чем надо понимать, что
большинство дорог будут
представлены линиями,

11.

Подготовка данных OSM к
map_polygons$geometry
Последний этап
анализу
=
подготовки это
map_polygons$geometry
%>% s2::s2_rebuild()
%>%
sf::st_as_sfc()
map_lines$geometry =
map_lines$geometry %>%
s2::s2_rebuild() %>%
sf::st_as_sfc()
map_polygons =
map_polygons %>%
st_collection_extract(
type="POLYGON")
map_polygons$natural[!
is.na(map_polygons$bui
lding)] = "yards"
map_polygons$natural
=
11
2022-03-15
as.factor(map_polygons
исправление геометрий
объектов загруженной
нами карты, который
необходим из-за
небольшого различия в
формате описания
полигонов в формате osm
и объектах класса sf.
Решается эта проблема
фукнциями
s2_rebuild(),st_as_sfc(), для
полигонов лучше также
добавить
st_collection_extract(type=“
POLYGON”). Для того
чтобы на картах получить
подписи на русском языке,

12.

Создание
собственных
карт
ggplot()+
geom_sf(map_polygons
%>%
filter(is.na(building)
),
map=aes(fill=natural))
+
geom_sf(map_polygons
%>% filter(!
is.na(building)),
map=aes(fill=building)
)+
geom_sf(map_lines,
map=aes(color=highway)
)+12
2022-03-15
theme_bw()

13.

Создание собственных
интерактивных
карт
c
leaflet() %>%
помощью
Leaflet
addProviderTiles("Open
TopoMap") %>%
addPolygons(data =
park_sf$geometry %>%
s2::s2_rebuild() %>%
sf::st_as_sfc(),
color = "black") %>%
addPolygons(data =
map_polygons,
fillColor = "grey",
color = "grey") %>%
addPolylines(data =
map_lines, color =
"red")
13
2022-03-15

14.

Создание собственных карт с
ggmap(mad_map)+
подложкой
geom_sf(data =
map_polygons,
map=aes(fill=natural),
inherit.aes = FALSE )+
geom_sf(data =
map_lines,
map=aes(color=highway)
, inherit.aes = FALSE)
+
xlab("Широта, °")+
ylab("Долгота, °")+
guides(fill =
guide_legend(title="Ле
генда"),
color2022-03-15
14
=guide_legend(title =

15.

Извлечение
данных
OSM
area =
После отображения
st_area(park_sf)
wood_area =
sum(st_area(map_polygo
ns %>%
filter(natural=="Древе
сные насаждения")))
water_area =
st_area(map_polygons
%>%
filter(natural=="Водны
е объекты"))
building_area =
sum(st_area(map_polygo
ns %>% filter(!
is.na(building))))
footway_length =
st_length(map_lines)
15
2022-03-15
footway_area =
данных, попробуем их
проанализировать и
составить собственую
сводную таблицу. Для
этого посчитаем ощую
площадь парка - st_area(),
площадь древесных
насаждений, водных
объектов и зданий. Для
объектов отображенных
линиями, мы можем
посчитать только их
длинну - st_length(), но
предположив, что ширина
всех тропинок
приблизительно 2 метра,
можем оценить их
площадь

16.

Рассчет
запечатанности
summary_ha = tibble(
Для создания свободной
water =
round(as.numeric(water
_area )/10000, 2),
wood =
round(as.double(wood_a
rea)/10000, 2),
build =
round(as.double(buildi
ng_area)/10000, 2),
road =
round(as.double(footwa
y_area)/10000, 2),
grass =
round(as.numeric(area)
/10000, 2) - water wood - build -road,
area =
16
2022-03-15
round(as.numeric(area)
таблицы переведем
данные в гектары,
посчитаем долю площади
каждого типа объекта от
общей площади и
рассчитаем степень
запечатанности, как
отношений площади
территории под зданиями
и дорожками к общей
площади, без учета
водных объектов.

17.

Создание
сводной
таблицы
library(gt)
С помощью пакета gt
#library(gtsummary)
library(flextable)
library(huxtable)
summary_table =
rbind(summary_ha,summa
ry_perc,
summary_final)
summary_table =
summary_table %>%
select(name,water,road
,wood,build,grass,area
)
sum_tabl =
summary_table %>%
gt(rowname_col =
"name")
%>%
17
2022-03-15
tab_header(
соберем результаты в
сводную таблицу

18.

Итог
18
2022-03-15