Наука о данных

И. В. Щуров, НИУ ВШЭ

На странице курса находятся дополнительные материалы.

Домашнее задание №14

Автор задания Анна Денисова.

За домашнее задание всего можно набрать 25 баллов, хотя для получения полного балла достаточно набрать 20 баллов. Однако, в отличие от предыдущих заданий, баллы, набранные сверх 20, также будут учтены (то есть можно набрать «больше максимума», и в итоге получить больше 10 баллов за курс).

Сейчас мы решим практическую задачу, касающуюся работы с геоданными, чтобы было понятно, где это вообще может пригодиться. Для начала скачаем данные для задания. Перейдите по ссылке и скачайте архив, который содержит полигоны (так мы будем называть многоугольники в координатах земной поверхности) административно-территориальных единиц по всей России в формате GeoJSON. Из данного архива вам понадобится файл admin_level_8.geojson.

Давайте посмотрим, что он в себе содержит:

In [ ]:
import json

with open('admin_level_8.geojson', encoding = 'utf-8') as f:
    a = json.load(f)
    
type(a)
In [ ]:
a.keys()

Как понятно из названия, формат GeoJSON представляет геоданные в виде JSON-файла. Сами полигоны представлены списком, который хранится по ключу "features". Давайте посмотрим на первый элемент этого списка:

In [ ]:
a['features'][0]

Мы видим, что данный элемент имеет тип Multipolygon, то есть какое-то объединение полигонов, однако фактически это обычный полигон. Как так получилось? По сути, это список полигонов, состоящий из одного элемента. Разработчики выбранной базы данных скорее всего сделали это, чтобы привести данные к единому формату. Преобразуем этот список координат в объект Polygon из библиотеки shapely.

In [ ]:
from shapely.geometry import Polygon

poly = Polygon(a['features'][0]['geometry']['coordinates'][0][0]) 
## Мы берем первый элемент из списка полигонов
poly

Давайте также на всякий случай проверим, что внутри нашего файла только объекты типа Multipolygon:

In [ ]:
arr = []
for i in range(len(a['features'])):
    arr.append(a['features'][i]['geometry']['type'])

print(set(arr))

(Задание 1: 0.5 балла) У каждого полигона есть центр. Получите координаты центра полигона выше с помощью метода .centroid. Выведите широту и долготу численно.

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 2: 0.5 балла) Для наших последующих целей нам понадобятся только те полигоны, центроиды которых лежат в центральных районах Москвы. Импортируйте файл Полигон Москвы.geojson так же, как мы импортировали admin_level_8.geojson. В этом файле лежит единственный полигон (список features имеет длину 1). Вытащите из GeoJSON-файла координаты единственного лежащего в нем полигона, преобразуйте в Polygon и сохраните в переменную moscow_poly. Обратите внимание, что теперь вам придется работать не с форматом Multipolygon, а с форматом Polygon, поэтому вам не нужно будет опускаться в рамках списка координат на уровень ниже (нужно брать не a['features'][0]['geometry']['coordinates'][0][0], а a['features'][0]['geometry']['coordinates'][0]).

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 3: 2 балла) Создайте датафрейм, колонками которого будут id объектов из a, их названия, а также широта и долгота центра мультиполигонов из файла admin_level_8.geojson. Назовите датафрейм df, а колонки с широтой и долготой — lat и lon соответственно. Это важно для последующих заданий (задание 5).

Hint: для создания датафрейма необходимо в цикле пройтись по массиву features из a.

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 4: 0.5 балла) Создайте датафрейм из файла с полигоном Москвы. В нем должна быть одна строка и два столбца: в одном название ("Москва"), в другом - полигон moscow_poly. Назовите датафрейм df1, а колонку с полигоном poly - это важно для последующих заданий (задание 5).

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 5: 3 балла) Мы хотим понять, какие из центров полигонов датафрейма df лежат внутри полигона Москвы. Это поможет нам выделить полигоны центральных районов Москвы. Для этого мы воспользуемся функцией библиотеки GeoPandas sjoin. Заметим, что sjoin работает только для геодатафреймов, поэтому сначала нам нужно создать два геодатафрейма из df и df1. Посмотрим, как это делается:

In [ ]:
import geopandas as gpd

gdf = gpd.GeoDataFrame(df, geometry = gpd.points_from_xy(df['lon'], df['lat']))
# при создании геодатафрейма нужно явно указывать, какой столбец содержит информацию
# о геоданных (geometry)
# в случае с gdf мы сказали, что этот столбец нужно создать, используя
# координаты точек, хранящиеся в столбцах `lon` и `lat`

gdf1 = gpd.GeoDataFrame(df1, geometry = 'poly')
# в этом случае мы сказали, что нужно просто взять содержимое столбца `poly`
# в качестве `geometry`
In [ ]:
gdf.head(2)
In [ ]:
gdf1.head(2)

Примените функцию sjoin с параметром op="intersects" к геодатафреймам выше и создайте геодатафрейм только с теми полигонами, центроиды которых лежат в полигоне Москвы. Прочитайте подробнее о spatial join в документации. Обратите внимание на параметр how и вспомните пары по SQL. Что вам нужно выбрать: right/left/inner?

Hint: специфицировать никакие параметры, кроме op и how, не нужно.

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 6: 2 балла) Выше вы видели пример создания геодатафрейма (задание 5). Создайте датафрейм с колонками id (id полигона), name (название полигона) и poly (геометрическая сущность Polygon, полученная с помощью библиотеки shapely). Преобразуйте данные в геодатафрейм. Возьмите данные из файла admin_level_8.geojson. В геодатафрейме должны быть только полигоны, центроиды которых лежат в полигоне Москвы, их id вы должны были получить в предыдущем задании. Обратите внимание, что теперь геометрией должен выступать полигон (параметр geometry в GeoDataFrame).

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 7: 1.5 балла) Соберите данные о пунктах выдачи заказов отсюда. Список точек загружается с помощью javascript-кода отдельным запросом, и вы можете перехватить ответ на этот запрос. Подробнее о том, как это сделать, было рассказано в этом видео. Вам нужно найти правильный адрес запроса (подсказка: он заканчивается буквой s). Затем просто обратиться к нему с помощью requests.get() и преобразовать результат в питоновский объект с помощью метода .json() (как мы это делали в ДЗ про API).

Hint \#1: Вы можете получить данные по колонкам id, coordinates, dtype, isWb,pickupType, longitude, latitude, status

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 8: 1 балл) Разбейте столбец с координатами на два столбца с широтой и долготой.

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 9: 3 балла) Создайте геодатафрейм, где геометрией выступают координаты точек WildBerries. Примените функцию sjoin к данному датафрейму и датафрейму, полученному в задании 6. Удалите из рассмотрения все точки WildBerries, координаты которых лежат не в центральных районах Москвы.

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 10: 3 балла) Посмотрите, в каких из центральных районов больше пунктов выдачи. Создайте геодатафрейм, в котором геометрией являются полигоны, и раскрасьте районы в зависимости от количества пунктов выдачи.

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿

(Задание 11: 6 баллов) Данное задание взято из домашки по машинному обучению ФКНа и адаптировано под наши цели.

Скачайте файл с данными об остановках общественного транспорта в формате .xlsx по ссылке (так гарантированно не возникнет проблем с парсингом файла) и загрузите таблицу в ноутбук. Если вдруг сайт Правительства Москвы сойдет с ума, то возьмите какую-нибудь версию данных отсюда. Для удобства визуализации мы будем работать только с остановками в ЦАО.

In [ ]:
data = pd.read_excel('City surface public transport stops.xlsx')
data = data[data.AdmArea_en == "Czentral`ny'j administrativny'j okrug"]
data = data.reset_index()
data.head()

Воспользуемся библиотекой folium для визуализации данных:

In [ ]:
import folium

m = folium.Map([55.75215, 37.61819], zoom_start=12)
for ind, row in data.iterrows():
    folium.Circle([row.Latitude_WGS84_en, row.Longitude_WGS84_en],
                  radius=10).add_to(m)
m

Как вы уже могли заметить, для каждой остановки указаны номера маршрутов, проходящих через неё. Логично соединить ребрами соседние остановки каждого маршрута. Однако мы не знаем, в каком порядке автобусы объезжают остановки. Но мы можем применить эвристический алгоритм, который восстановит нам порядок маршрутов:

  • Для каждого маршрута выделим список всех остановок, через которые он проходит.
  • Выберем начальную остановку маршрута как точку, наиболее удаленную от всех остальных остановок этого маршрута.
  • Каждую следующую точку маршрута будем выбирать как самую близкую из оставшихся точек (не включенных в маршрут ранее).

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

In [ ]:
def get_routes(data):
    '''
    Accumulate routes from raw data
    param data: pd.DataFrame - public transport stops data
    return: dict - unsorted stops ids for each route,
                   e.g. routes['A1'] = [356, 641, 190]
    '''

    # YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿
    raise NotImplementedError


def sort_routes(data, routes):
    '''
    Sort routes according to the proposed algorithm
    param data: pd.DataFrame - public transport stops data
    param routes: dict - unsorted stops ids for each route
    return: dict - sorted stops ids for each route
    '''

    # YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿
    raise NotImplementedError
In [ ]:
routes = get_routes(data)
sorted_routes = sort_routes(data, routes)

(Задание 12: 2 балла) Визуализируйте любые 10 из предложенных маршрутов на карте с помощью библиотеки folium. Для этого воспользуйтесь методом folium.vector_layers.PolyLine().

In [ ]:
# YOUR CODE HERE ‿︵‿︵ヽ(°□° )ノ︵‿︵‿