Привет, Хабр!

На связи участник профессионального сообщества NTA Максим Алёшин.

Сегодня создание интерактивной карты на python не составляет большого труда: стоит подключить библиотеку (например, Folium или Bokeh), указать картографический сервер, и после выполнения нескольких «магических» строк кода ваши данные как на ладони!

Но что делать, если данные есть, визуализировать их хочется, а в сети, в которой вы работаете, нет доступа к картографическим серверам? В этом случае помогут Plotly и GeoPandas. Но придётся провести ряд подготовительных работ. В этом посте шаг за шагом я покажу, как построить интерактивную карту России по регионам с помощью Plotly, которая будет работать без интернета, регистрации и смс.

Введение

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

Plotly — мощная библиотека визуализации данных на Python, которая позволяет создавать широкий спектр интерактивных визуализаций, включая карты. Одним из преимуществ Plotly является то, что она работает с объектами других библиотек Python, таких как Pandas и NumPy. В Plotly «из коробки» есть классы для работы с картами (см. примеры), но они также требуют доступа к внешним серверам, мне же такое не подходит по условию задачи.

Помимо стандартных модулей типа Pandas и NumPy будет необходима GeoPandas. Это библиотека с открытым исходным кодом, упрощающая работу с геопространственными данными в Python за счёт расширения типов данных, используемых в Pandas. Она позволяет выполнять пространственные операции над геометрическими типами. Другими словами, это те же панды, но знающие толк в картографии. GeoPandas оперирует с геометрическими объектами, работа с которыми реализована в библиотеке Shapely.

Прежде чем начать

Первым делом я реализую простой класс для создания интерактивной карты России. Это будет своего рода слой‑подложка, который заменит существующий в Plotly класс plotly.graph_objects.layout.Geo, работающий с подготовленным слоем карты.

Разберу, почему не подходит готовый инструмент:

• подложка карты грузится с внешнего сервера;

• набор возможных для отображения областей: «africa», «asia», «europe», «north america», «south america», «usa», «world»;

• отсутствие «красивых» систем координат для отображения России.

Первый недостаток для моей задачи ­‑ главный. Из второго следует, что с координатами России мне пришлось бы работать на карте всего мира. Здесь имеются отдельные слои для Азии, Европы, но наша широкая и необъятная страна распростёрлась на обе части света. Также по умолчанию в Plotly карты рисуются в стандартной цилиндрической системе координат (широта и долгота). Есть и другие проекции, но они далеки от идеальных для двумерной проекции нашей страны на плоскости. Когда же нужно визуализировать данные непосредственно на карте России, хочется, чтобы она выглядела как в контурных картах.

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

Основные чекпоинты реализации:

• подготовить картографические данные для основного слоя карты;

• реализовать класс для слоя-подложки карты;

• построить интерактивную карту по регионам, изобразить города на карте;

• добавить на карту информацию в виде всплывающих подсказок.

Весь код с реализованным модулем для построения карты России и используемые данные доступны в репозитории.

Описание данных

Чтобы не отвлекаться на подготовку данных для последующей визуализации, в репозитории приведены уже подготовленные датасеты. Географические данные представлены в формате GeoJSON.

GeoJSON — формат для кодирования различных структур географических данных, базирующийся на известном формате JSON. Он поддерживает такие типы геометрии, как: Point, LineString, Polygon, MultiPoint, MultiLineString и MultiPolygon.

В моём арсенале имеются два предобработанных датасета:

• С границами регионов (russia_regions.geojson): географические данные, актуальные на 2021 год, взяты из репозитория. Имеются поля с названием региона, федерального округа, добавлено поле с населением. Данные о населении извлечены со страницы Население субъектов Российской Федерации Википедии.

• С населением городов (russia_cities_population.parquet): имеются поля с названием города, региона, федерального округа, численностью населения и географическими координатами (взяты из стороннего датасета).

Предобработка географических данных

Процесс предобработки данных для создания слоя подложки карты реализован в Jupyter ноутбуке geodata_preparation.ipynb.

Импорт библиотек

Для предобработки географических данных импортирую уже упомянутые библиотеки и matplotlib для пробной визуализации.

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from tqdm import tqdm

import geopandas as gpd
import shapely
from shapely.ops import snap, unary_union
from shapely.geometry.polygon import Polygon
from shapely.geometry.multipolygon import MultiPolygon

Считываю географические данные:

gdf = gpd.read_file("data/russia_regions.geojson")
gdf.info()            
gdf.sample(4)

Результат:

<class 'geopandas.geodataframe.GeoDataFrame'>
RangeIndex: 85 entries, 0 to 84
Data columns (total 4 columns):
 #   Column            Non-Null Count  Dtype   
---  ------            --------------  -----   
 0   region            85 non-null     object  
 1   federal_district  85 non-null     object  
 2   population        85 non-null     int64   
 3   geometry          85 non-null     geometry
dtypes: geometry(1), int64(1), object(2)
memory usage: 2.8+ KB

Каждая строка здесь содержит информацию об одном регионе. Последний столбец — самый главный. Собственно, именно наличием столбца geometry GeoPandas датафреймы отличаются от пандасовских. В случае с границами объектов значения этого столбца имеют тип данных Polygon или MultiPolygon. В первом случае объект можно описать одним простым многоугольником. Второй тип описывает элемент, состоящий из нескольких многоугольников. Он хранится в виде массива Polygon. Большинство регионов России имеют в составе острова или по другим причинам состоят из нескольких разрозненных территорий.

Системы координат

С помощью GeoPandas можно с лёгкостью построить геометрию. Посмотрим, как будет выглядеть карта в том виде, как она есть в датасете.

gdf.plot()
plt.title('Default. Cylindric CRS');

Первая проблема, которая бросается в глаза — пересечение со 180 меридианом откидывает часть Чукотки. Вторая — искажение геометрии из‑за цилиндрической системы координат. Хочется видеть нашу родную постройней и пособранней как‑то. Для этого достаточно перейти в другую картографическую систему координат. И это легко сделать с помощью GeoPandas.

Разработчики GeoPandas добавили в библиотеку датасет геодезических параметров EPSG, в котором содержатся параметры перехода из одной географической системы отсчёта в другую и связанные с ними единицы измерения. Каждый код этого имеет значение в диапазоне от 1024 до 32 767. Также сюда добавили и параметры некоторых систем отсчёта ESRI. Стандартная привычная система отсчёта в единицах широты и долготы имеет код EPSG:4326.

Я перепробовал довольно много систем отсчёта и в итоге нашёл ту, которая сделает Россию снова красивой. Проверить существующие коды EPSG или ESRI можно на epsg.io. Ниже привожу выборку из 4 систем.

crss = ['EPSG:3576', 'EPSG:5940', 'ESRI:102027', 'EPSG:32646']
fig = plt.figure()
for i, crs in enumerate(crss):
    ax = fig.add_subplot(2, 2, i+1)
    gdf.to_crs(crs).plot(ax=ax)
    ax.set_title(crs)
plt.tight_layout()
plt.show()

Последняя прямо как из контурных карт! Её и оставлю.

gdf = gdf.to_crs('EPSG:32646')

По горячим следам нужно разобраться с раздробленностью Чукотки в датасете. Для этого пробегу по циклу по всем парам полигонов и объединю их с допуском в 100 единиц с помощью команды .snap(). Построю геометрию до и после преобразования.

CHUK = 'Чукотский автономный округ'

gdf.loc[gdf.region == CHUK].plot(figsize=(6,6), facecolor="none")
plt.show()

# Объединение разбитых полигонов Чукотки
new_chuk = []
chuk_geoms = gdf.loc[gdf.region == CHUK, 'geometry'].values[0].geoms

# приклеим друг к другу полигоны находящиеся на расстоянии менее 100 единиц
for i, g in enumerate(chuk_geoms):
    new_g = g
    for j in range(len(chuk_geoms)):
        new_g = snap(new_g, chuk_geoms[j], 100)
    new_chuk.append(new_g)
new_chuk = unary_union(MultiPolygon(new_chuk))
gdf.loc[gdf.region == CHUK, 'geometry'] = new_chuk

gdf.loc[gdf.region == CHUK].plot(figsize=(6,6), facecolor="none")
plt.show();

Упрощение геометрии

Одна из основных проблем построения интерактивных карт – размерность геометрии. Чем «детальней» объекты, тем дольше будет идти прорисовка. После нескольких попыток изобразить границы регионов как есть на интерактивном виджете Plotly, стало ясно, что точек многовато. Было решено: упрощать! Ведь взаимодействие с картой становится намного приятней, если она не «лагает». Да и в последующем нужно будет наносить на карту новые объекты.

Основным инструментом оптимизации будет выступать специальный метод .simplify(), который позволяет уменьшить количество точек. Он принимает на вход параметр tolerance, все точки упрощённой геометрии будут находиться на расстоянии не более этого допуска. Рассмотрю уменьшение на примере Калининградской области. Упрощу геометрию региона до расстояния между соседними точками в 500 единиц.

i = 61 # Калининградская область
tol = 500
print(f"Точек {shapely.get_num_coordinates(gdf.geometry[i])}") 
display(gdf.geometry[i])
print(f'Точек {shapely.get_num_coordinates(gdf.geometry[i].simplify(tol))}') 
display(gdf.geometry[i].simplify(tol))

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

Ниже функция упрощения геометрии всех регионов. Использую модуль tqdm для отслеживания прогресса выполнения.

def prepare_regions(gdf, area_thr=100e6, simplify_tol=500):
    """Подготовка регионов к построению
    
    - Упрощение геометрии с допуском simplify_tol
    - Удаление полигонов с площадью менее area_thr
    """
    gdf_ = gdf.copy()
    
    # Вспомогательный столбец для упорядочивания регионов по площади
    gdf_['area'] = gdf_.geometry.apply(lambda x: x.area)

    # Удаляем маленькие полигоны
    tqdm.pandas(desc='Удаление мелких полигонов')
    gdf_.geometry = gdf_.geometry.progress_apply(lambda geometry: 
        MultiPolygon([p for p in geometry.geoms if p.area > area_thr]) 
            if type(geometry) == MultiPolygon else geometry
    )
    
    # Упрощение геометрии
    gdf_.geometry = gdf_.geometry.simplify(simplify_tol)
    
    geoms = gdf_.geometry.values
    pbar = tqdm(enumerate(geoms), total=len(geoms))
    pbar.set_description_str('Объединение границ после упрощения')
    # проходим по всем граничащим полигонам и объединяем границы
    for i, g in pbar:
        g1 = g
        for g2 in geoms:
            if g1.distance(g2) < 100:
                g1 = snap(g1, g2, 800)
        geoms[i] = g1
    gdf_.geometry = geoms
    
    # сортировка по площади
    gdf_ = gdf_.sort_values(by='area', ascending=False).reset_index(drop=True) 
    
    return gdf_.drop(columns=['area'])

prepare_regions выполняет сразу несколько задач, которые описаны несколькими блоками кода. Функция принимает датасет gdf, порог площади area_thr и значение допуска для упрощения геометрии simplify_tol. Полигоны со значением площади меньше area_thr будут удалены:

• Вначале создаётся вспомогательный столбец для упорядочивания регионов по площади, чтобы самые маленькие регионы рисовались последними. Например, если Адыгею рисовать раньше Краснодарского края, то она оказывается им закрашена.

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

• Упрощается геометрия. После этого необходимо объединить границы, чтобы на карте было как можно меньше «дыр».

Финальная обработка геометрии

Последнее, что необходимо сделать — преобразовать полигоны и мультиполигоны шэйпы, которые способна считывать для построения Plotly. Графические объекты поддерживают рисование GeoJSON объектов только на упомянутых слоях типа go.layout.Geo, от которых я отказался. Поэтому остаётся только использовать класс go.Scatter, а значит, необходимо перевести геометрии в формат отдельных массивов координат x и y. Для этого я реализовал функцию‑преобразование geom2shape.

def geom2shape(g):
    """Преобразование полигонов и мультиполигонов в plotly-readable шэйпы    
    
    Получает на вход Polygon или MultiPolygon из geopandas, 
    возвращает pd.Series с координатами x и y
    """
    # Если мультиполигон, то преобразуем каждый полигон отдельно, разделяя их None'ами
    if type(g) == MultiPolygon:
        x, y = np.array([[], []])
        for poly in g.geoms:
            x_, y_ = poly.exterior.coords.xy
            x, y = (np.append(x, x_), np.append(y, y_))
            x, y = (np.append(x, None), np.append(y, None))
        x, y = x[:-1], y[:-1]
    # Если полигон, то просто извлекаем координаты
    elif type(g) == Polygon:      
        x, y = np.array(g.exterior.coords.xy)
    # Если что-то другое, то возвращаем пустые массивы
    else:
        x, y = np.array([[], []])
    return pd.Series([x,y])

Если преобразуется мультиполигональный объект, то последовательность координат каждого полигона отделяется от последующего с помощью None. Так plotly в режиме заливки будет закрашивать каждый замкнутый контур отдельно, но при этом на карте это будет считаться единым объектом. Не баг, а фича! Осталось применить написанные функции к датасету, и сохранить его.

# Упрощение геометрии   
regions = prepare_regions(gdf)
# Преобразование полигонов в шейпы   
regions[['x','y']] = regions.geometry.progress_apply(geom2shape)
# Запись на диск   
regions.to_parquet('data/russia_regions.parquet')

Слой-подложка карты

Всё готово, чтобы создать слой‑подложку для карты. Построение буду проводить с помощью собственного класса mapFigure, который наследуется от класса фигуры Plotly — go.Figure. Для этого импортирую подмодуль plotly.graph_objects, который содержит основные классы для построения фигур. В листинге ниже приведён код модуля map_figure.py. Помимо класса mapFigure в нём реализована функция convert_crs для преобразования координат в двух массивах x_arr, y_arr из одной системы отсчёта в другую. Она пригодится для подготовки данных к отображению на карте. Замечу, что у go.Figure разработчики отключили возможность добавления новых атрибутов и методов. Но здесь можно обойтись переопределением метода инициализации, а дальше работать с картой, как с обычным go.Figure.

'''Класс для слоя подложки карты России'''

import pandas as pd
import geopandas as gpd
import plotly.graph_objects as go
from shapely.geometry import Point

REGIONS = pd.read_parquet("data/russia_regions.parquet")

def convert_crs(x_arr, y_arr, to_crs='EPSG:32646', from_crs="EPSG:4326"):
    """Преобразование значений координат в массивах x_arr и y_arr
    из географической системы отсчёта from_crs в систему to_crs
    """
    data = [Point(x,y) for x,y in zip(x_arr, y_arr)]
    pts = gpd.GeoSeries(data, from_crs).to_crs(to_crs)
    
    return pts.x, pts.y

class mapFigure(go.Figure):
    """ Шаблон фигуры для рисования поверх карты России
    """
    def __init__(self, # дефолтные параметры plotly
        data=None, layout=None, frames=None, skip_invalid=False, 
        **kwargs # аргументы (см. документацию к plotly.graph_objects.Figure())
    ):
        # создаём plotlу фигуру с дефолтными параметрами
        super().__init__(data, layout, frames, skip_invalid, **kwargs)

        # прорисовка регионов
        for i, r in REGIONS.iterrows():
            self.add_trace(go.Scatter(x=r.x, y=r.y,
                                      name=r.region,
                                      text=r.region,
                                      hoverinfo="text",
                                      line_color='grey',
                                      fill='toself',
                                      line_width=1,
                                      fillcolor='lightblue',
                                      showlegend=False
            ))
        
        # не отображать оси, уравнять масштаб по осям
        self.update_xaxes(visible=False)
        self.update_yaxes(visible=False, scaleanchor="x", scaleratio=1)

        # чтобы покрасивее вписывалась карта на поверхности фигуры
        self.update_layout(showlegend=False, dragmode='pan',
                           width=800, height=450, 
                           margin={'l': 10, 'b': 10, 't': 10, 'r': 10})

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

• x,y — подготовленные ранее шэйпы для построения;

• name — название объекта;

• text — информация, отображаемая при наведении на объект;

• hoverinfo ‑ определяет, какая информация будет отображаться при наведении. В моём случае будет отображаться информация, прописанная в параметре text;

• line_color, line_width — цвет и ширина границ;

• fill — определяет, как именно объекты будут заполнятся цветом; в режиме «toself» регионы окрашиваются полностью до своих границ;

• fillcolor — устанавливает цвета заливки;

• showlegend — определяет, показывать ли легенду.

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

Работа с объектом mapFigure

Примеры построения карты можете найти в Jupyter ноутбуке russia_map.ipynb.

Сразу импортирую написанный модуль и сделаю другие необходимые импорты.

from map_figure import mapFigure, convert_crs

import plotly.graph_objects as go
import plotly.express as px
import pandas as pd

Вид по умолчанию

Создаю объект карты и отображаю её.

russia_map = mapFigure()
russia_map.show()

Результат: стандартный объект Plotly go.Figure c возможностью перемещаться по фигуре, зуммировать. Из интерактива пока только отображение названия региона. При реализации класса я намеренно использовал единственное поле датасета. Так моделирую ситуацию, когда в моём распоряжении только названия географических областей и их контуры. Поскольку названия этих областей хранятся в параметре name каждого объекта фигуры, то есть возможность изменять их.

Изменение вида основного слоя

В примере ниже поменяю цвет каждого региона, отображу его название и численность населения в нём, название федерального округа во всплывающей подсказке. Цвет определю уникальным для каждого федерального округа. Для этого использую пресеты палитр из px.colors.qualitative.

import pandas as pd
import plotly.express as px

regions = pd.read_parquet("data/russia_regions.parquet")
fo_list = list(regions['federal_district'].unique())
colors = px.colors.qualitative.Pastel1

for i, r in regions.iterrows():
    popul_text = f"Население: <b>{r.population:_} </b>".replace('_', ' ')
    text = f'<b>{r.region}</b><br>{r.federal_district} ФО<br>{popul_text}'
    russia_map.update_traces(selector=dict(name=r.region),
        text=text,
        fillcolor=colors[fo_list.index(r.federal_district)])
russia_map.show()

Добавление элементов

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

df = pd.read_parquet("data/russia_cities_population.parquet")
df.sample(4)

Не забываю перевести координаты городов в нужную систему отсчёта с помощью convert_crs. Теперь добавляю точки на фигуру также, как я делал бы это в go.Figure.

# координаты в требуемой для карты системе отсчёта
df['x'], df['y'] = convert_crs(df.lon, df.lat) 
# палитра цветов для точек
cities_palette = {fd: px.colors.qualitative.Dark2_r[i] for i, fd in enumerate(df.federal_district.unique())}

# добавление точек на фигуру
russia_map.add_trace(go.Scatter(
    x=df.x, y=df.y, name='города',
    text="<b>"+df.city+"</b><br>Население <b>"+df.population.astype('str')+"</b>",
    hoverinfo="text", showlegend=False, mode='markers',
    marker_size=df.population/1e4, marker_sizemode='area', marker_sizemin=3,
    marker_color=df.federal_district.map(cities_palette)
))

russia_map.show()

Приближу:

Как видно дальнейшая работа ничем не отличается от стандартной работы с Plotly.

Заключение

Итак, интерактивная карта, работающая без доступа к картографическим серверам, — не миф, а реальность. Замечу, что помимо границ регионов можно визуализировать и природные географические объекты (реки, озёра) города, населённые пункты и прочее. В свободном доступе есть картографические данные с детализацией до уровня улиц, в том числе в GeoJSON формате. Да, в созданном здесь шаблоне детализация невысока, но его можно улучшать. Здесь я продемонстрировал подход и показал, что, несмотря на ограничения сети, с помощью Plotly можно достичь впечатляющих результатов.

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

Спасибо за прочтение! Буду рад комментариям и идеям по улучшению.