knitr::opts_chunk$set(echo = TRUE)
#в основном здесь нужен будет dplyr (для работы с дата фреймом) и ggplot для графики
library(tidyverse)
#пакет для чтения таблиц excel
library(readxl)
#пакет для того, чтобы загрузить данные с http://www.naturalearthdata.com/
#загрузка в формате shape файлов
#нужен будет, чтобы загрузить границы стран
library(rnaturalearth)

#пакет для работы с пространственными данными
#сейчас используется намного меньше, чем следующий пакет, но все же используется, потому что от него зависят некоторые другие пакеты
library(sp)

#основной и самый используемый сейчас пакет для работы с векторными пространственными данными
library(sf)

#пакет для создания тематических карт
library(tmap)

#сборник цветовых палитр
library(RColorBrewer)

Первая часть (карты с помощью пакета ggplot2): https://rpubs.com/baltti/thematic-maps-ggplot

Вторая часть (карты с помошью пакета cartography): https://rpubs.com/baltti/thematic-maps-cartography

Данные будут те же, что и в предыдущих частях.

Сначала подготовка данных. В файле, скачанном с http://happyplanetindex.org/ лучше поправить заголовки: часть заголовков идет в 2 абзаца в исходном документе, а также убрать пробелы в заголовке (из-за этого могут выдаваться ошибки)

countries<-ne_countries(scale = "medium", continent = "europe",
                        type = "countries", returnclass = "sf") %>% 
  st_crop(c(xmin = -10, xmax = 45, 
            ymin = 34, ymax = 71))%>%
  st_transform(3857)
## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries
hpi<- sp::merge(countries, 
                read_excel("hpi-data-2016.xlsx",col_names = TRUE,
                     sheet = 5,range = "B6:O146"), 
        by.x="sovereignt",by.y="Country")%>% 
  select(1,64:77)

Тематические карты с помощью пакета tmap

Общая информация о пакете, а также ссылки на ресурсы по использованию пакета: https://github.com/mtennekes/tmap

Вводная статья по использованию пакета: https://cran.r-project.org/web/packages/tmap/vignettes/tmap-getstarted.html

И есть статья авторов пакета: https://www.jstatsoft.org/article/view/v084i06

Кроме того в книге Geocomputation in R глава в основном построена на использовании этого пакета.

Пакет предназначен специально для создания тематических карт: tmap - thematic maps. Отрисовка также как и в пакете cartography производится слоями.

Особенность пакета в собственном внутреннем синтаксисе, отчасти похожем на синтаксис пакета ggplot2. Оба эти пакета основаны на приципах grammar of graphics.

Но главная особенность и преимущество этого пакета, на мой взгляд, в том что он позволяет создавать интерактивные карты из статических заменой одной строчки кода. Интерактивные карты здесь создаются на базе js-библиотеки leaflet.

Вообще для R есть одноименный пакет leaflet, но с ним работать не всегда удобно, потому что в нем есть некоторые особенности добавления объектов (часто требуется широта и долгота, не всегда sf объекты корректно наносятся).

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

Основные функции

В пакете есть основные функции отрисовки (main plotting method), у каждой из которых свои настройки (aesthetics):

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

Между собой слои могут быть скомбинированы, чтобы добавить следующий слой к предыдущему, нужно просто поставить + между строками (аналогично ggplot2). Но важно помнить, что перед любой из этих функций должна идти функция tm_shape, аналогично тому, как в пакете cartography сначала нужно добавить слой с геометрией объектов. Это нужно, чтобы задать проекцию, границы слоев (bounding box) и, конечно, же исходную геометрию.

Дополнительно у этим еще есть функции настройки дополнительных элементов.

Кроме перечисленных есть еще функция создания “быстрых” тематических карт - qtm.

С помощью пакета tmap можно делать как статические карты, так и интерактивные. Тип карты задается функцией tmap_mode.

Картограмма или choropleth map

В этом пакете три варианта отрисовки полигонов: полигон с внутренней заливкой и границами, только внутренняя заливка и только границы полигона.

Создадим самую простую карту с полигонами стран.

tm_shape(hpi)+
  tm_polygons()

# в этом случае можно воспользоваться tm_sf, результат будет аналогичный

Далее мы можем настраивать эти полигоны так как нам нужно. Зададим цвет полигонов в зависимости от величины индекса счастья.

tm_shape(hpi)+
  tm_polygons(col="Happy_Planet_Index", #переменная в зависимости от величины которой будет задаваться цвет
              style="equal", n=4, #метод классификации и количество классов
              palette="Pastel2")+
  tm_layout(legend.outside = TRUE, #легенда будет размещена рядом с картой, а не внутри
            legend.outside.position = "right") 

Аналогичная карта, но только с внутренней заливкой полигонов.

Классификацию здесь можно задать как с использованием встроенных опций классификации

tm_shape(hpi)+
  tm_fill(col="Happy_Planet_Index", #переменная в зависимости от величины которой будет задаваться цвет
              style="equal", n=4, #метод классификации и количество классов
              palette="Pastel2") +
  tm_layout(legend.outside = TRUE, #легенда будет размещена рядом с картой, а не внутри
            legend.outside.position = "right") 

Также классификацию можно задавать вручную

tm_shape(hpi)+
  tm_polygons(col="Happy_Planet_Index", 
              style="fixed", #при таком стиле классификации интервалы задаются вручную
              breaks = c(10,20,30,40), 
              palette="Pastel2",title="Индекс счастья",
              lwd=1, lty="dashed") + #толщина и тип линии для границ
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "right")

Эти карты статические, но их легко превратить в интерактивные.

tmap_mode('view') #интерактивная карта
## tmap mode set to interactive viewing
tm_shape(hpi)+
  tm_fill(col="Happy_Planet_Index", 
              style="jenks", n=4, 
              palette="Pastel2")+
  tm_borders(col="#2966a3",lwd = 2)

Для интерактивной карты можно задавать параметры, например, подложки.

tmap_mode('view')
## tmap mode set to interactive viewing
tm_shape(hpi)+
  tm_polygons(col="Happy_Planet_Index", style="jenks", n=4, palette="RdYlGn",
              title="Индекс счастья",
              alpha = 0.4) + #прозрачность
  tmap_options(basemaps = 'Stamen.TonerLite')

Пропорциональные символы или bubble map

Для создания карты пропорциональных символов можно воспользоваться несколькими вариантами функций. Основная функция здесь tm_symbols. Кроме нее еще есть tm_bubbles, tm_dots, tm_markers, tm_squares - фактически они отличаются только самим типом символа, аргументы для всех этих функций аналогичны.

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

#извлечение центроидов
centroid<-st_centroid(hpi, of_largest_polygon = TRUE)
## Warning in st_centroid.sf(hpi, of_largest_polygon = TRUE): st_centroid assumes
## attributes are constant over geometries of x
tmap_mode('plot') #статическая карта
tm_shape(hpi)+
  tm_borders(lwd = 1)+
tm_shape(centroid)+
  tm_symbols(size="Happy_Planet_Index")

Настроим символы и добавим подписи к ним.

centroid<-st_centroid(hpi, of_largest_polygon = TRUE)
## Warning in st_centroid.sf(hpi, of_largest_polygon = TRUE): st_centroid assumes
## attributes are constant over geometries of x
tmap_mode('plot')
## tmap mode set to plotting
tm_shape(hpi)+
  tm_borders(lwd = 1)+
tm_shape(centroid)+
  tm_symbols(size="GDP",col = "#167b9c",
             scale=1.2,border.col = NA,
             title.size = "ВВП, $/чел.")+ 
  tm_text(text = "sovereignt",size = 0.7,just = "top")+
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "right")

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

tmap_mode('plot')
## tmap mode set to plotting
tm_shape(hpi)+
  tm_borders(col="#167b9c", lwd = 1)+
tm_shape(centroid)+
  tm_symbols(shape= "Happy_Planet_Index", shapes.style = "jenks",
             shapes.n=4, shapes = c(pch = 21:25), 
             col = "#45a16b",title.shape = "Индекс счастья") + 
  tm_text(text = "sovereignt",size = 0.5,just = "left")+
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "right")

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

tmap_mode('view')
## tmap mode set to interactive viewing
tm_shape(centroid)+
  tm_bubbles(size= "Happy_Planet_Index", 
             col="GDP", style="jenks", n=5, pal="RdYlGn", 
             alpha = 0.7) +
  tmap_options(basemaps = 'Esri.WorldTerrain')
## Legend for symbol sizes not available in view mode.

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

При этом, если карта интерактивная, то есть возможность эти слои включать и отключать (под кнопками увеличения/уменьшения слева)

tmap_mode('view')
## tmap mode set to interactive viewing
tm_shape(hpi)+
  tm_polygons(col="Happy_Planet_Index", style="jenks", n=4, palette="RdYlGn",
              title="Индекс счастья", alpha = 0.4) +
tm_shape(centroid)+
  tm_bubbles(size="GDP",col = "#5229a3", alpha = 0.4)
## Legend for symbol sizes not available in view mode.

“Быстрые” тематические карты

Этот вариант позволяет создавать карты минимум кода одной функцией. Но здесь намного меньше возможностей настройки вида карты

tmap_mode('plot')
## tmap mode set to plotting
qtm(hpi, fill = "Happy_Planet_Index", style = "col_blind")+
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "right")

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

tmap_mode('plot')
## tmap mode set to plotting
qtm(hpi, fill = "Happy_Planet_Index", style = "col_blind",by="Region")+
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "right")

Компоновка нескольких карт

В пакете существуют встроенные функции компоновки tmap_arrange и tm_facets. Первая функция позволяет компоновать созданные карты на одну картинку.

equal<-tm_shape(hpi)+
  tm_polygons(col="Happy_Planet_Index", style="equal", n=4, palette="Pastel2",
              title="Индекс счастья")+ #здесь title это заголовок легенды
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "bottom", main.title = "Равные интервалы",
            main.title.size = 1)

quant<-tm_shape(hpi)+
  tm_polygons(col="Happy_Planet_Index", style="quantile", n=4, palette="Pastel2",
              title="Индекс счастья")+
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "bottom", main.title = "Квантили",
            main.title.size = 1)

jenks<-tm_shape(hpi)+
  tm_polygons(col="Happy_Planet_Index", style="jenks", n=4, palette="Pastel2",
              title="Индекс счастья")+
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "bottom", main.title = "Естественные интервалы",
            main.title.size = 1)

tmap_arrange(equal,quant,jenks) 

Функция tm_facets предназначена для компоновки нескольких карт двумя способами: 1. группировка объектов по одной или двум переменным (например, по названиям континентов или стран, или по величинам какой-то переменной) 2. создание нескольких карт на основе нескольких переменных

tmap_mode('plot')
## tmap mode set to plotting
tm_shape(hpi)+
  tm_polygons()+
  tm_facets(by="Region", #переменная по которой будет проводиться группировка объектов
            free.coords = FALSE) #должны ли для каждой карты быть свои пределы координат

Это работает и с интерактивными картами

tmap_mode('view')
## tmap mode set to interactive viewing
tm_shape(hpi)+
  tm_polygons(col = "Happy_Planet_Index", style="jenks", n=4, palette="Pastel2",
              title="Индекс счастья")+
  tm_facets(by="Region", #переменная по которой будет проводиться группировка объектов
            free.coords = FALSE, #должны ли для каждой карты быть свои пределы координат
            sync = FALSE) #должны ли интерактиные панели быть синхронизированы

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

tmap_mode('view')
## tmap mode set to interactive viewing
tm_shape(hpi)+
  tm_polygons(col = "Happy_Planet_Index", style="jenks", n=4, palette="Pastel2",
              title="Индекс счастья")+
  tm_facets(by="Region", as.layers = TRUE) 

Также можно компоновать карты по двум переменным для этого в аргументе by нужно указать вектор из необходимых переменных. Кроме того с помощью аргумента along можно разделить объекты, чтобы отрисовывать карты на нескольких картинках.

tmap_mode('plot')
## tmap mode set to plotting
tm_shape(hpi)+
  tm_polygons(col = "#66aacc")+
  tm_facets(by="sovereignt", along = "Region",
            drop.empty.facets =  TRUE,
            free.coords = TRUE) 

Из полученных нескольких карт можно создать анимацию встроенной в пакет функцией tmap_animation (требует установки ImageMagick)

Сравнение пакетов

Для сравнения создадим одну и ту же картограмму с помощью трех разных пакетов: ggplot2, cartography и tmap.

#ggplot2
library(classInt) #для разбивки на интервалы
int<-classIntervals(hpi$Happy_Planet_Index,n=4, style = "jenks")
hpi$interval <- cut(hpi$Happy_Planet_Index, 
                    breaks = int$brks,
                    include.lowest = FALSE)
hpi %>%
  ggplot() + 
  geom_sf(aes(fill=interval))+ 
  scale_color_brewer(type = "div",palette="RdYlGn",  
                     aesthetics = "fill", name = "Индекс счастья", 
                     guide = guide_legend(keyheight = unit(5, units = "mm"),
                     title.position = 'top', reverse = T))+  
  theme_minimal()+
  ggtitle("ggplot2 package")

#cartography
library(cartography)
par(mar = c(0, 0, 1.2, 0))
plot(hpi$geometry, col = "white", border = "black")
choroLayer(hpi, var = "Happy_Planet_Index", 
           col = brewer.pal(n = 4,name = "RdYlGn"), 
           method = "fisher-jenks", nclass = 4, 
           border = "black", lwd = 0.5,legend.pos = "right", 
           legend.title.txt = "Индекс счастья", add = T)
layoutLayer(title = "Cartography package",  frame = FALSE, scale = NULL, 
            theme = "green.pal", north = FALSE)

#tmap
tmap_mode('plot')
## tmap mode set to plotting
tm_shape(hpi)+
  tm_polygons(col="Happy_Planet_Index", style="jenks", n=4, palette="RdYlGn",
              title="Индекс счастья")+
  tm_layout(legend.outside = TRUE, 
            legend.outside.position = "right", main.title = "tmap package",
            main.title.size = 1)

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

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

У пакета tmap, к сожалению, не очень много вариантов создания тематических карт, только самые основные, но его большое преимущество - это возможность создавать интерактивные карты минимумом усилий. Да, конечно, есть и другие пакеты для этого, например, mapview, plotly и leaflet, но они проигрывают в простоте настройки визуализации.