Делаем хитмап стоимости недвижимости на общедоступных алгоритмах и опенсорсных библиотеках

uedufy_xbmeezfo2imfanm84dmk.png

Карта 2gis.ru работает на WebGL-движке, который позволяет визуализировать данные. Когда мы делали слой недвижимости, то решили добавить ещё и тепловую карту стоимости квадратного метра. Во-первых, это красиво. А во-вторых, такие карты могут помочь с выбором квартиры.

Под катом — про то, какими алгоритмами пользовались и с какими трудностями встретились, когда делали хитмап стоимости квартир на собственных данных.
В самом начале мы задались классическими вопросами:

  1. Какие вообще бывают тепловые карты?
  2. Есть ли готовые решения и подойдут ли они нам?
  3. Решает ли какой-нибудь сервис подобную задачу и как?

Вот что мы выяснили.

Типы хитмапов


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

Тепловые карты плотности — карты, где цвета кодируют плотность данных. Подходят для визуализации количественных значений: землетрясений, ДТП, плотности населения и т. д.

image
Карта частоты землетрясений на базе Mapbox

Тепловые карты значений — карты, где цвета кодируют абсолютное значение параметра в точке. С их помощью можно отобразить качественные значения какой-нибудь величины в зависимости от координат — например, температуру воздуха или скорость ветра.

image
Тепловая карта температуры на сервисе windy.com

Готовые решения


Для хитмапов плотности существует довольно много различных JS-библиотек, в том числе и с использованием WebGL. Можно выделить Simpleheat, Heatmap.js, WebGL Heatmap, WebGL Mapbox Heatmap, которые довольно шустро рассчитываются на клиенте по исходному набору данных.

Есть ещё heatmapper, который позволяет генерировать векторные хитмапы плотности. Правда, он написан на R, что нам не очень подходит. Вот репозиторий.

Перечисленные библиотеки хитмапов плотности основаны на алгоритме ядерной оценки плотности.

С хитмапами значений все немного сложнее. По сути, гуглится лишь одна JS-библиотека Temperature Map, реализующая интерполяционный алгоритм Inverse distance weighting (IDW). И ещё есть вариант реализации IDW на R.

В отличие от карт плотности, расчёт хитмапов значений занимает достаточно много времени. Поэтому их обычно генерируют заранее и отдают на клиент уже в качестве готовых растровый изображений или векторных геометрий.

Эксперименты


Использовать хитмапы плотности для отображения абсолютных значений не получится, потому что результат зависит от количества данных, а не от значения в точках. Этот вывод вполне очевиден — его можно сделать сразу после прочтения о типах тепловых карт. И можно было бы приступить к реализации хитмапа значений. Но вы даже не представляете, сколько раз меня спрашивали, а почему бы не заиспользовать хитмап плотности, вот у Мапбокса он прям на лету рисуется. В итоге решили показать разницу наглядно.

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

Для экспериментов взяли данные о средней стоимости продажи квадратного метра жилой недвижимости в зданиях Москвы. Их рассчитали по данным объявлений наших партнёров.

Строим тепловую карту плотности. Написали код в RStudio на основе реализации в heatmapper. Стоимость квадратного метра указали в качестве веса для функции kde2d.weighted. Получилась такая картина.

image
Хитмап плотности на данных о стоимости продажи недвижимости

По задумке хитмап должен характеризовать стоимость квадратного метра: зеленый — дёшево, желтый — средне, красный — дорого.

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

image
Искусственно уплотнённые данные на хитмапе плотности породили новое «дорогое» пятно

Мы ожидаем, что на хитмапе не будет принципиальных изменений, ведь в среднем цена в этом квадрате не изменилась, но тёмно-красное пятно говорит об обратном.

Теперь пробуем хитмап значений. В RStudio набросали код генерации с использованием алгоритма IDW. Результат на тех же данных получился более интересным и логичным.

image
Хитмап значений на данных о стоимости продажи недвижимости

Уплотнили данные — и тепловая карта практически не изменилась. По крайней мере, не появилось новых цветовых пятен.

rog0leyh31t65v-lkkpgiynl7wi.png
Слева исходные данные, а справа — искусственно уплотнённые так же,
как и для хитмапа плотности

Изменения, конечно, присутствуют. Но они связаны с тем, что из-за увеличенного количества данных изменились интервалы стоимости, характеризующие цвет — потому что интервалы не фиксированы, а рассчитываются путём кластеризации всех значений, присутствующих в данных.

Хитмапы других сервисов


На просторах интернета нашли две отличных реализации хитмапов стоимости недвижимости — Яндекс.Недвижимость и Квартиры-Домики.рф.

941vcvydij_0nddivns4ifxeefy.png
Тепловая карта стоимости недвижимости на Яндекс.Недвижимость

_vndr5kmttgszskxvaims3mlzui.png
Тепловая карта стоимости продажи недвижимости на «Квартирах-Домиках»

Про алгоритмы расчёта хитмапов у Яндекса нашлась только статья, что они его сделали, и немного про то, как подготавливались данные для алгоритма интерполяции. А вот про реализацию тепловой карты на «Квартирах-Домиках» есть полезная техническая статья.

Нам понравились оба варианта. В решении Яндекса отметили для себя интересный внешний контур хитмапа — по сути, это контур плотности данных. Он выглядит вполне логичным, ведь независимо от выбора интерполяционного алгоритма точность предсказаний будет падать с удалением от точек с известными значениями. В решении на «Квартирах-Домиках» приглянулись округлые контуры ценовых зон, потому что такие контуры выглядят более естественными.

Реализация


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

Первым делом сформулировали требования к внешнему виду тепловых карт.

Требования


  1. Шкала для тепловых карт должна быть дискретной, то есть сам хитмап должен иметь ограниченное количество цветов и каждый цвет должен характеризовать некоторый интервал значений.
  2. Контуры цветовых зон должны быть округлыми, так как это выглядит естественнее.
  3. Тепловая карта должна быть обрезана по контуру плотности данных, потому что чем дальше мы находимся от известных значений, тем точность интерполяции хуже.

Из экспериментов на R стало ясно, что два первых пункта закроет алгоритм IDW вместе с кластеризацией множества известных значений. А контур плотности по исходному набору данных построить вообще не проблема — с этим, как минимум, поможет D3. И если с D3 всё понятно — берём да используем, — то с IDW стоит разобраться.

Что такое IDW


Суть алгоритма IDW такова: значения в произвольной точке рассчитываются на основе взвешенных значений в точках исходных данных. Вес зависит от расстояния до интерполируемой точки, чем расстояние меньше, тем вес больше, и наоборот. Таким образом, точки исходных данных, расположенные рядом с интерполируемой точкой, будут больше влиять на итоговое значение, чем точки, расположенные далеко.

Базовая формула IDW выглядит так:

$ u(x)={\begin{cases}{\dfrac {\sum _{i=1}^{N}{w_{i}(x)u_{i}}}{\sum _{i=1}^{N}{w_{i}(x)}}}&{\text{если }}d(x, x_{i})\neq 0{\text{ для всех }}i,\\u_{i}&{\text{если }}d(x ,x_{i})=0{\text{ по крайней мере для одного }}i,\end{cases}}$


где:

Увеличение степени алгоритма $p$ будет повышать влияние точек данных, расположенных рядом, и уменьшать влияние удалённых. При больших значениях получится что-то похожее на диаграмму Вороного.

Есть ещё вариант модифицированного вычисления веса, при котором для расчёта значения в точке используются только данные, находящиеся в радиусе $R$ от интерполируемой точки. Выглядит эта формула так:

$w_{k}(x)=\left({\frac {\max(0, R-d(x, x_{k}))}{Rd(x, x_{k})}}\right)^{2}$

Вместе с k-мерным деревом или подобной ему структурой, которая реализует быстрый пространственный поиск, эта модификация даёт серьёзный прирост в производительности.

Для стоимости недвижимости вариант с радиусом выглядит более подходящим, ведь стоимость жилья в отдалённых районах вряд ли существенно влияет на стоимость в интерполируемой точке. Плюс такой вариант просто быстрее. Поэтому на нём и остановились.

Собираем алгоритм


Результатом работы IDW будет двумерный массив значений в каждой точке области, по которой у нас есть данные — например, в городе. Сопоставив каждому значению цвет из палитры, мы получим растровое изображение хитмапа, которое можно наложить на карту.

Но карту можно зумить, а вот растровые изображения нельзя — при масштабировании они сильно теряют в качестве. Поэтому мы рассчитывали хитмап по тайлам в соответствии с тайловой сеткой интересующей нас области.

А ещё нужно обрезать тепловую карту по контуру плотности. Контур нужно строить на всём массиве исходных данных, а сам хитмап рассчитывается по тайлам. Значит, нам придется порезать контур плотности на тайлы и учитывать соответствующие тайлы контура плотности при расчёте тайлов тепловой карты. С этим нам поможет библиотека geojson-vt.

В итоге пришли к такому алгоритму:

  1. Рассчитываем баунд (минимальный прямоугольник, содержащий все точки) исходных данных.
  2. С каждой стороны добавляем к нему отступы размером, равным радиусу R для алгоритма IDW, тем самым расширяя баунд.
  3. Считаем тайловую сетку для расширенного баунда.
  4. Строим контур плотности данных и режем его на тайлы в соответствии с тайловой сеткой баунда.
  5. Рассчитываем значение пикселей каждого тайла в баунде — используем алгоритм IDW с модифицированной формулой расчёта весов с радиусом R. Значения пикселей вне соответствующего контура плотности оставляем пустыми.
  6. Для быстрого поиска точек в радиусе R от рассчитываемого пикселя используем эту реализацию k-мерного дерева.

Что получилось


Реализуем алгоритм на нашем любимом JavaScript и прогоняем на данных о стоимости продажи квадратного метра жилой недвижимости в Москве.

0ipvwmrkpr3yydv4eh8t4nrful8.png

Хитмап получился нормальным. Но выглядел так, будто заболел ветрянкой: на нём было много очень мелких областей-пятен. Чтобы понять, в чём дело, мы отобразили точки из набора данных, по которым рассчитывался хитмап.

y3rzhxme-wdh5hqve4yqg0nmlyw.png
Участок хитмапа с отображёнными точками из набора данных

Оказалось, что мелкие области формируются вокруг близко расположенных точек со значениями из разных цветовых интервалов. Чтобы избавиться от этого, надо как-то уменьшить плотность данных. Для этого построили равномерную квадратную сетку по баунду данных и рассчитали медианные значения для каждого квадрата из тех, в которые попали хотя бы n точек. Остальные квадраты просто проигнорировали.

Такой результат получился для квадратной сетки со стороной квадрата 500 м и n = 2.

wf9njbnditdug3ahhgza1zasec8.png

Этот вариант нам понравился больше. На нём и остановились.

Приятный бонус — время расчёта хитмапа после группировки данных значительно уменьшилось. На одной и той же выборке данных расчёт тайлов со второго по 12-й зум для Москвы занимал 28 минут на сырых данных и 2 минуты на сгруппированных.

Так к алгоритму построения хитмапов в самое начало добавился ещё один пункт — предварительно сгруппировать данные по квадратной сетке.

Результат


Мы построили тепловые карты продажи и аренды жилой недвижимости в десяти городах России. Они доступны в онлайн-версии 2ГИС в режиме поиска недвижимости. Чтобы было удобнее посмотреть и сравнить, собрали список:

Планы


Наша карта — трёхмерная и векторная. Её можно вращать, наклонять и непрерывно зумить. А хитмапы пока что растровые и сгенерированы для дискретных целочисленных зумов. Получается, например, что на зуме 13.9 мы используем тайлы хитмапа с 13-го зума, увеличивая их пропорционально. Из-за этого контуры цветовых зон выглядят недостаточно гладкими, несмотря на двойную плотность тайлов.

4ek98nm43dncg7hzmyp1e14-gig.png

Чтобы решить эту проблему, планируем сделать тайлы хитмапов векторными. Для этого будем генерировать на основе растровых тайлов векторные контуры и дорабатывать движок карты.

P. S.


Мы ни в коем случае не претендуем на то, что изобрели тепловые карты. Статья лишь про то, как с помощью общедоступных алгоритмов и опенсорсных библиотек мы визуализировали свои данные о продаже и аренде недвижимости на нашей карте.

Всем дешёвых квадратных метров :)

© Habrahabr.ru