Объединение открытых данных Open Street Map и Landsat для уточнения площадей зеленых зон

Моя цель - предложение широкого ассортимента товаров и услуг на постоянно высоком качестве обслуживания по самым выгодным ценам.

Привет, Хабр! Под катом мы хотели бы поговорить о том какие преимущества привносит в нашу жизнь объединение пространственных данных из различных открытых источников. Рассматривать, для примера, будем следующую задачу: оценить насколько в “зеленом” районе расположен тот или иной объект недвижимости. Без долгих предисловий - переходим к делу!

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

Подход 1. Самый точный, но не для ленивых

Представим, что такая задача по оценке возникла для объекта города Санкт-Петербург по адресу “27 к1, Комендантский проспект, округ Юнтолово, Санкт-Петербург, Северо-Западный федеральный округ, 197371, Россия”. Задачу получил не какой-то там штатный программист, а, предположим, человек не умеющий программировать, но старательный и готовый изучать что-то новое. Поискав что тут можно придумать в интернете, “специалист по оценке зелености района” уже скачивает QGIS и учится пользоваться плагином Quick Map Services. Следующие шаги могут быть очевидны, но отнюдь не просты: нужно разобраться с тем в каком виде представляются пространственные данные в геоинформационных системах (ГИС) (векторный и растровый), затем следует научиться создавать такие объекты и сопоставлять друг с другом (чтобы площади посчитать).

Рисунок 1. Картинка чтобы не было скучно читать: вот такие примитивы используются в ГИС, чтобы описать пространственные объекты
Рисунок 1. Картинка чтобы не было скучно читать: вот такие примитивы используются в ГИС, чтобы описать пространственные объекты

Несмотря на то что попутно с описанными выше действиями придется восполнять множество других пробелов в знаниях (например о том что такое система координат и какие они бывают), самым времязатратным процессом может оказаться этап создания объектов - векторизация (когда буквально обводишь те объекты, которые видно на снимке). И если со знаниями и опытом операции в ГИС будут каждый раз даваться все легче и легче, то вот векторизация тяжело поддается автоматизации (как правило).

Итак, что же может сделать наш специалист - он может загрузить данные спутника Google Satellite через модуль Quick Map Services и вручную разметить данные - самостоятельно нарисовать геометрии объектов. В данном случае - парков и все что на них похоже. Затем сможет для простоты объяснения подхода наложить буфер в 2000 метров вокруг объекта недвижимости посчитать его площадь. Затем выделить те объекты, которые лежат в том же районе и выделены им как зеленые зоны, и рассчитать их площадь. Отношение “площадь зеленых зон в районе / площадь района” и будет искомым значением.

Важное уточнение - предложим, что наш специалист живет в Санкт-Петербурге и прекрасно знает исследуемый район. Поэтому вручную выделенные объекты по снимку можно будет назвать некоторой комбинацией наземной рекогносцировки вместе с ручным дешифрированием.

Анимация 1. Векторизация по простому. Выделение объектов на карте и представление их в виде векторного слоя
Анимация 1. Векторизация по простому. Выделение объектов на карте и представление их в виде векторного слоя

Преимущества: Достаточно точная и объективная оценка. Рассчитанное отношение в полной мере характеризует “зеленость района” - наш специалист определенно получит премию

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

Подход 2. Дешево и сердито

Ну конечно же речь про OpenStreetMap… А пока представим, что наш лирический герой понял, что IT в целом намного интереснее сферы недвижимости и начал изучать язык программирования Python (а какой же ещё). Он познакомился с модулем osmnx и библиотекой geopandas для обработки пространственных данных и библиотекой shapely для обработки геометрий. При помощи этой троицы уже можно процедурно проводить оценку зелености территории. Для этого нужно: 

  1. Создать полигон вокруг точки - границы района для анализа

  2. Запросить данные OSM, которые лежат в этом полигоне (для запроса следует разобраться с системой тегов, которая есть в OSM)

  3. Посчитать площади у полигона района и полученных геометрий парков и зеленых зон OSM

Такой подход значительно быстрее, так как не требует ручного создания полигонов. Действительно, за нас уже эти полигоны создали другие люди - этим и прекрасен OSM. Однако OSM имеет и свои недостатки - данные могут быть неточны. Более того, даже если полигоны пользователями отрисованы корректно, всегда присутствует вероятность, что вы просто прогадаете с набором тегов для запроса и пропустите какой-нибудь важный кусок данных (Рисунок 2).

Рисунок 2. Пример зеленой зоны, которая по данным OSM является дворовой территорией и по тегам хоть как-то относящимся к “зеленым зонам” выделить её никак не получается 
Рисунок 2. Пример зеленой зоны, которая по данным OSM является дворовой территорией и по тегам хоть как-то относящимся к “зеленым зонам” выделить её никак не получается 

И таких зон довольно много. На Рисунке 2 показана лишь небольшая часть. Однако уже такой подход позволяет быстро и, хотя примерно, но всё же объективно оценить зеленость территории.

Преимущества: Быстрый и простой способ оценки. Несмотря на неточности, данные OSM в большинстве случаев использовать можно

Недостатки: Высокая погрешность в оценке, если некорректны исходные данные, или если набор тегов по которым осуществляется запрос, составлен недостаточно подробно

Подход 3. Снимки Landsat как источник объективных данных

О, сколько всего было сказано и написано про Landsat. А так же про приложения на основе рассчитанного вегетационного индекса NDVI. Если заинтересовались, то предлагаем изучить следующие материалы: 

  • NDVI - теория и практика

  • Вегетационные индексы

  • Как фермеру узнать состояние своих полей по NDVI?

  • В целом про продукты Landsat: Диапазоны Landsat 8 в работе

Если же высказаться коротко, то: 

  1. Landsat - это такой спутник, который летает над Землей и делает снимки с довольно высоким пространственным разрешением (от 15-30 в оптическом спектре до 100 метров в дальнем инфракрасном)

  2. Normalized Difference Vegetation Index (NDVI) - это вегетационный индекс, который можно рассчитать на основе снимка, сделанного Landsat. NDVI принимает значения от -1 до 1, где -1 означает, что в данном пикселе совсем нет ничего зеленого (например, вода), а 1 значит что цветут луга и колышутся березы. 

И да, при использовании NDVI стоит помнить про его недостатки. Довольно подробно этот вопрос разбирался на одном из семинаров спбгеотеха - см. запись “Эдуард Казаков - Блеск и нищета NDVI (#спбгеотех)”.

Итак, мы можем получить снимок Landsat (скачать “руками” можно тут - ищите продукт “Landsat 8-9 OLI/TIRS C2 L2”) для выбранной территории и рассчитать NDVI. Затем, варьируя порог мы можем сегментировать изображение на два класса: 0 - не зеленая зона и 1 - зеленая зона. Мы можем варьировать порог как хотим (Рисунок 3), но, честно говоря, наверняка не знаем какой порог бинаризации подойдет лучше всего. Более того, такой “оптимальный порог” будет отличаться для каждого нового снимка полученного с Landsat: как для новой территории, так и за отличный от текущего момент времени.  

Рисунок 3. Сегментация матрицы NDVI по порогам 0.15, 0.2, 0.25. Черным цветом показаны границы “зеленых зон”
Рисунок 3. Сегментация матрицы NDVI по порогам 0.15, 0.2, 0.25. Черным цветом показаны границы “зеленых зон”

Преимущества: Подход позволяет получать актуальную информацию о том насколько территория зеленая

Недостатки: Выбор оптимального значения порога затруднен - не понятно какое значение выбрать, чтобы по порогу отделить растительность не от растительности. Есть конечно источники, в которых можно найти информацию о предлагаемых порогах, например от 0.2 и выше. Но стоит держать в голове, что это оптимальное значение порога может изменяться от снимка к снимку, от сезона к сезону и т.д.

Подход 4. Площади по Landsat, пороги - по OpenStreetMap

Логичным решением проблемы с выбором порога может стать использование данных OSM. Действительно, как мы заметили раньше, даже несмотря на то, что в данных OSM могут быть артефакты или отсутствовать часть реально существующих объектов, в целом они дают представление о том где какие объекты (парки, скверы, здания) находятся. Так давайте же наложим на карто-схеме геометрии полученные по OSM и данные NDVI, а потом построим гистограмму (Рисунок 4). 

Рисунок 4. Частотная гистограмма по двум классам: значения NDVI в погонах “зеленых объектов” по данным OSM и остальные значения NDVI на снимке
Рисунок 4. Частотная гистограмма по двум классам: значения NDVI в погонах “зеленых объектов” по данным OSM и остальные значения NDVI на снимке

На рисунке черной линией показано значение, которое можно использовать как порог. Тогда можно дополнить данные OSM дополнительными полигонами, полученными при помощи автоматической векторизации матрицы Landsat. Таким образом, все что расположено справа от рассчитанного порога превращается в “зеленые зоны”. За счет этого происходит расширение границ (Рисунок 5).

Рисунок 5. Границы зеленых зон до и после уточнения по данным Landsat
Рисунок 5. Границы зеленых зон до и после уточнения по данным Landsat

Вместо заключения

Будет картинка…

Рисунок 6. Сравнение выделенных границ объектов различными способами
Рисунок 6. Сравнение выделенных границ объектов различными способами

Как видно, третий подход с уточнением границ объектов по данным Landsat оказался наиболее приближенным к эталону (ручной векторизации). Цифры площади тоже это подтверждают: 

  • (Эталон) Ручная векторизация: 28.3%

  • Расчет по данным OSM: 13.9%

  • Расчет по уточненным данным OSM и NDVI Landsat: 34.9% 

Как видим, наша оценка хотя и немного завышена, она все равно намного ближе к действительности, чем расчет по данным OSM, который в два раза меньше. 

P.S. На основе данной идеи мы оставили следующие артефакты: 

  • https://github.com/red5ai/estaty - Python open-source библиотека для обработки пространственных данных и подготовки прототипов алгоритмов для анализа объектов недвижимости. Библиотека совсем молодая, но мы планируем её развивать и улучшать

  • http://api.greendex.wiredhut.com/ - сервис для рассчета “зеленого” индекса и внутренний инструмент (расширение для google spreadsheets), чтобы можно было удобно работать с вышеописанными алгоритмами как через API, так и через формулы  в таблицах (для тех, кто не умеет программировать). Данная версия использует упрощенный расчет без использования спутников Landsat

С рассказом про объединения данных OSM и Landsat был Михаил Сарафанов и Wiredhut team

Источник: https://habr.com/ru/articles/764686/


Интересные статьи

Интересные статьи

Всем привет! Я Станислав Бушуев, Software Engineer в Semrush. Сегодня хочу поделиться идеями, как можно реализовать синхронизацию данных между различными хранилищами. Такие задачи иногда возникают в р...
Мир онлайн-покупок становится всё привычнее, а значит, и обезличенных данных про каждого пользователя всё больше. Билайн ТВ использует для онлайн-кинотеатра рекомендательную систему на основе данных: ...
Всем привет! Cегодня я расскажу, как применять InfluxDB для мониторинга систем хранения данных. Я затрону следующие темы: Наше путешествие в стек мониторинга InfluxDB для мониторинга систем хр...
Если на вашем сайте присутствует большое количество контента, то для отображения пользователю его приходится так или иначе делить.Все известные мне способы имеют недостат...
Привет, Хабр! Сегодня хотим представить вам некоммерческий открытый датасет, собранный командой студентов магистратуры «Наука о данных» НИТУ МИСиС и Zavtra.Online (подразделении SkillFact...